Merge pull request #2 from mowcom/equation-update

.cursor/rules/agent.md

@@ -0,0 +1,5 @@
+The agent should only request explicit user permission to complete a suggested code change if the change is large, potentially destructive, or could have significant side effects. For most routine, incremental, or clearly safe changes, the agent should proceed and apply the change without asking for additional user approval.
+
+Summary of intent:
+- Small, routine, or obviously safe changes: Apply directly, no need to ask.
+- Large, sweeping, or potentially risky changes: Ask for user confirmation before proceeding. 

.cursor/rules/er_model_guidelines.md

@@ -28,7 +28,7 @@ globs:
     - black
     - isort
   ```
-
+y
 # Repository Structure
 - `src/` → core modules (e.g. `er_model.py`, `allometry.py`, `metrics.py`).
 - `scripts/` → CLI entrypoints (e.g. `run_pipeline.py`, `analyze_results.py`).

.cursorignore

@@ -0,0 +1,5 @@
+# .cursorignore
+node_modules/
+outputs/
+*.csv
+*.ipynb

.gradio/certificate.pem

@@ -0,0 +1,31 @@
+-----BEGIN CERTIFICATE-----
+MIIFazCCA1OgAwIBAgIRAIIQz7DSQONZRGPgu2OCiwAwDQYJKoZIhvcNAQELBQAw
+TzELMAkGA1UEBhMCVVMxKTAnBgNVBAoTIEludGVybmV0IFNlY3VyaXR5IFJlc2Vh
+cmNoIEdyb3VwMRUwEwYDVQQDEwxJU1JHIFJvb3QgWDEwHhcNMTUwNjA0MTEwNDM4
+WhcNMzUwNjA0MTEwNDM4WjBPMQswCQYDVQQGEwJVUzEpMCcGA1UEChMgSW50ZXJu
+ZXQgU2VjdXJpdHkgUmVzZWFyY2ggR3JvdXAxFTATBgNVBAMTDElTUkcgUm9vdCBY
+MTCCAiIwDQYJKoZIhvcNAQEBBQADggIPADCCAgoCggIBAK3oJHP0FDfzm54rVygc
+h77ct984kIxuPOZXoHj3dcKi/vVqbvYATyjb3miGbESTtrFj/RQSa78f0uoxmyF+
+0TM8ukj13Xnfs7j/EvEhmkvBioZxaUpmZmyPfjxwv60pIgbz5MDmgK7iS4+3mX6U
+A5/TR5d8mUgjU+g4rk8Kb4Mu0UlXjIB0ttov0DiNewNwIRt18jA8+o+u3dpjq+sW
+T8KOEUt+zwvo/7V3LvSye0rgTBIlDHCNAymg4VMk7BPZ7hm/ELNKjD+Jo2FR3qyH
+B5T0Y3HsLuJvW5iB4YlcNHlsdu87kGJ55tukmi8mxdAQ4Q7e2RCOFvu396j3x+UC
+B5iPNgiV5+I3lg02dZ77DnKxHZu8A/lJBdiB3QW0KtZB6awBdpUKD9jf1b0SHzUv
+KBds0pjBqAlkd25HN7rOrFleaJ1/ctaJxQZBKT5ZPt0m9STJEadao0xAH0ahmbWn
+OlFuhjuefXKnEgV4We0+UXgVCwOPjdAvBbI+e0ocS3MFEvzG6uBQE3xDk3SzynTn
+jh8BCNAw1FtxNrQHusEwMFxIt4I7mKZ9YIqioymCzLq9gwQbooMDQaHWBfEbwrbw
+qHyGO0aoSCqI3Haadr8faqU9GY/rOPNk3sgrDQoo//fb4hVC1CLQJ13hef4Y53CI
+rU7m2Ys6xt0nUW7/vGT1M0NPAgMBAAGjQjBAMA4GA1UdDwEB/wQEAwIBBjAPBgNV
+HRMBAf8EBTADAQH/MB0GA1UdDgQWBBR5tFnme7bl5AFzgAiIyBpY9umbbjANBgkq
+hkiG9w0BAQsFAAOCAgEAVR9YqbyyqFDQDLHYGmkgJykIrGF1XIpu+ILlaS/V9lZL
+ubhzEFnTIZd+50xx+7LSYK05qAvqFyFWhfFQDlnrzuBZ6brJFe+GnY+EgPbk6ZGQ
+3BebYhtF8GaV0nxvwuo77x/Py9auJ/GpsMiu/X1+mvoiBOv/2X/qkSsisRcOj/KK
+NFtY2PwByVS5uCbMiogziUwthDyC3+6WVwW6LLv3xLfHTjuCvjHIInNzktHCgKQ5
+ORAzI4JMPJ+GslWYHb4phowim57iaztXOoJwTdwJx4nLCgdNbOhdjsnvzqvHu7Ur
+TkXWStAmzOVyyghqpZXjFaH3pO3JLF+l+/+sKAIuvtd7u+Nxe5AW0wdeRlN8NwdC
+jNPElpzVmbUq4JUagEiuTDkHzsxHpFKVK7q4+63SM1N95R1NbdWhscdCb+ZAJzVc
+oyi3B43njTOQ5yOf+1CceWxG1bQVs5ZufpsMljq4Ui0/1lvh+wjChP4kqKOJ2qxq
+4RgqsahDYVvTH9w7jXbyLeiNdd8XM2w9U/t7y0Ff/9yi0GE44Za4rF2LN9d11TPA
+mRGunUHBcnWEvgJBQl9nJEiU0Zsnvgc/ubhPgXRR4Xq37Z0j4r7g1SgEEzwxA57d
+emyPxgcYxn/eR44/KJ4EBs+lVDR3veyJm+kXQ99b21/+jh5Xos1AnX5iItreGCc=
+-----END CERTIFICATE-----

README.md

@@ -1,153 +1,107 @@
-# ER Model - Nigeria Mangroves
+---
+title: Mangrove ER Model Dashboard
+emoji: 🌱
+colorFrom: green
+colorTo: blue
+sdk: gradio
+sdk_version: "5.29.0"
+app_file: dashboard/app.py
+pinned: false
+---
 
-A Python implementation of the Emissions Reduction (ER) model for Nigeria mangrove projects. This repository provides tools to estimate carbon sequestration over a 30-year period based on mangrove growth and survival rates.
+# Mangrove Emissions Reduction Model & Dashboard
 
-## Recent Updates and Improvements
-
-### Growth Model Enhancements (2024)
-- **Updated Allometric Equations**: Implemented Zanvo et al. (2023) species-specific equations:
-  - Rhizophora: Total = 1.938 × (DBH² H)^0.67628
-  - Avicennia: Total = 1.486 × (DBH² H)^0.55864
-
-- **Chapman-Richards Implementation**:
-  - Improved growth parameter estimation
-  - Direct use of growth rate parameter (b) instead of deriving from Tm
-  - Separate parameterization for DBH and height growth
-  - Current parameters:
-    * Rhizophora: max DBH 11.07cm, max height 12.0m
-    * Avicennia: max DBH 17.0cm, max height 8.8m
-
-- **Growth Rate Parameters**:
-  - b = 0.25 yr^-1 (independently estimated growth rate)
-  - c = 2.1 (shape parameter)
-  - Removed Tm-based calculation (previously b = ln(c)/Tm)
-
-### Dashboard Improvements
-- Added height growth visualization
-- Improved number formatting with comma separators
-- Enhanced scenario comparison features
-- Added biomass per tree tracking
-- Updated parameter input ranges
-
-### Configuration Updates
-- Restructured parameter files for clarity
-- Added detailed parameter documentation
-- Improved error handling for invalid inputs
-- Added support for scenario analysis
+This project implements a flexible, modular dashboard for mangrove carbon sequestration and emissions reduction, supporting multiple growth models. **Currently, only the Declining Increment model is active; all other models are archived for reference.**
 
 ## Features
+- Interactive Gradio dashboard for project simulation and analysis
+- Updated: grouped accordions, tooltips as labels, reset-to-defaults button, and compact two-column mortality layout
+- Modular growth model architecture (see `er_model_core/growth_models/`)
+- All dashboard tables and plots are driven by model outputs (no redundant calculations)
+- Editable planting schedule and all model parameters for Declining Increment model
+- Robust, transparent, and reproducible calculations
 
-- Species-specific growth modeling using Chapman-Richards equations
-- Biomass estimation using allometric equations
-- Carbon conversion and emissions reduction calculations
-- Interactive dashboard for parameter exploration using Gradio
-- Configurable planting schedules and mortality rates
-- Scenario comparison for different DBH ranges and growth rates
-- Per-species carbon tracking and cumulative hectare calculations
-
-## Model Parameters
-
-### Species Parameters
-- **Rhizophora**
-  - Initial values: DBH = 1.0cm, Height = 0.5m
-  - Chapman-Richards parameters: a = 11.07cm, b = 0.35, c = 2.1
-  - Configurable planting density and mortality rates
-
-- **Avicennia**
-  - Initial values: DBH = 1.0cm, Height = 0.5m
-  - Chapman-Richards parameters: a = 17.0cm, b = 0.35, c = 2.1
-  - Configurable planting density and mortality rates
-
-### Carbon Parameters
-- Biomass to carbon conversion factor
-- Carbon to CO2 conversion factor
-- Buffer pool percentage
-- Leakage percentage
-- Baseline emissions rate
-
-## Installation
+## Usage
 
-1. Install [Miniconda](https://docs.conda.io/en/latest/miniconda.html) or [Miniforge](https://github.com/conda-forge/miniforge)
+1. **Access the Dashboard (Recommended):**
+   - The Mangrove ER Model Dashboard is hosted online.
+   - **Access it here:** [https://huggingface.co/spaces/mowcom/er-model](https://huggingface.co/spaces/mowcom/er-model)
 
-2. Create and activate the environment:
+2. **(Optional) Run Locally for Development:**
    ```bash
+   # Using conda (recommended)
    conda env create -f environment.yml
    conda activate er-model
-   ```
-
-3. Install pre-commit hooks:
-   ```bash
-   pre-commit install
-   ```
-
-## Project Structure
-
-- `src/` - Core model implementation
-  - `er_model.py` - Main ER calculation logic
-  - `allometry.py` - Species-specific allometric equations
-  - `metrics.py` - Carbon conversion and statistics
-- `scripts/` - Command-line tools
-  - `run_pipeline.py` - Execute full ER calculation workflow
-  - `analyze_results.py` - Generate reports and visualizations
-- `dashboard/` - Gradio web application for interactive exploration
-- `configs/` - Default parameter files
-- `outputs/` - Generated results and figures
-- `tests/` - Test suite
-
-## Dashboard Features
-
-The interactive Gradio dashboard allows users to:
-1. Adjust species-specific parameters:
-   - Planting density (trees/ha)
-   - Year-by-year mortality rates
-2. Configure planting schedule:
-   - Area planted per year (ha)
-   - Up to 5-year planting schedule
-3. Modify carbon parameters:
-   - Buffer pool percentage
-4. Explore scenarios:
-   - Test different DBH ranges
-   - Adjust growth rate factors
-   - Compare results across 1000 ha scenarios
-5. View results:
-   - Carbon sequestration over time
-   - Milestone year comparisons
-   - Per-species carbon tracking
-   - Cumulative hectare metrics
-   - Scenario comparison tables
-
-## Usage
-
-1. Configure parameters in `configs/params.yaml`
-
-2. Run the full pipeline:
-   ```bash
-   python scripts/run_pipeline.py --config configs/params.yaml
-   ```
-
-3. Launch the dashboard:
-   ```bash
    python dashboard/app.py
    ```
-   This will start a Gradio server and open the dashboard in your default web browser.
-
-## Development
-
-- Code style is enforced using `black` and `isort`
-- Type hints are required for all public functions
-- Tests are written using `pytest`
+   - To stop the local server, press `Ctrl+C` in the terminal.
+
+## Growth Models
+- **Declining Increment (ACTIVE):**
+  - Annual increment declines linearly to zero; total size is the sum of non-negative increments (never decreases or goes negative).
+  - Both discrete and continuous versions are available (toggle via config).
+- **Archived Models:**
+  - **Chapman-Richards** (see `er_model_core/growth_models/chapman_richards.py`)
+  - **Simple Linear** (see `er_model_core/growth_models/linear.py`)
+  - **Linear Plateau** (see `er_model_core/growth_models/linear.py`)
+  - **Jimenez-Lugo** (see `models/jimenez_lugo/README.md`)
+
+## File Structure
+- `dashboard/app.py` — Gradio dashboard UI and logic (now model-driven, with grouped accordions and tooltips)
+- `er_model_core/er_model.py` — Core model logic, growth model selection, and all calculations
+- `er_model_core/growth_models/` — Modular growth model implementations (active and archived)
+- `er_model_core/types.py` — Shared dataclasses for model config
+- `er_model_core/config_loader.py` — Config loading and parsing
+- `configs/` — Model configuration YAMLs
+- `models/` — Reference model folders and documentation
+- `tests/` — Unit tests for model logic and outputs
+
+## Model-Driven Dashboard
+- All per-year, per-species metrics (surviving trees, DBH, height, biomass per tree, total biomass) are calculated and stored in the model (`model.species_metrics`).
+- Dashboard tables and plots read directly from model outputs for speed, consistency, and maintainability.
+- No redundant calculations in the dashboard.
+
+## Adding or Re-enabling Models
+1. **To add a new model:**
+   - Implement the growth function in `er_model_core/growth_models/`.
+   - Add a config YAML to `configs/`.
+   - Add a tab/entry in `dashboard/app.py` (if you want it visible in the UI).
+2. **To re-enable an archived model:**
+   - Uncomment the relevant import and config entry.
+   - Ensure the model is listed in `MODEL_CONFIGS` in `dashboard/app.py`.
+   - Add or update documentation in `models/` as needed.
+
+## Directory Structure
+```
+er-chapman-richards/
+├── README.md                # This file
+├── environment.yml          # Conda environment
+├── dashboard/               # Gradio dashboard app
+├── er_model_core/           # Core model code
+│   ├── er_model.py
+│   ├── growth_models/
+│   ├── types.py
+│   └── config_loader.py
+├── configs/                 # Model parameter YAMLs
+├── models/                  # Model-specific docs
+│   ├── original/README.md   # Archived
+│   ├── jimenez_lugo/README.md # Archived
+│   └── declining_increment/README.md # Active
+├── tests/                   # Pytest suite
+└── ...
+```
+
+## Development & Testing
+- Code style: [black](https://black.readthedocs.io/), [isort](https://pycqa.github.io/isort/)
+- Type hints required for all public functions
+- Run tests with:
+  ```bash
+  pytest
+  ```
 - Pre-commit hooks ensure code quality
 
-## Dependencies
-
-Key dependencies (specified in `environment.yml`):
-- numpy=1.24.3
-- pandas=2.0.3
-- gradio (latest version)
-- matplotlib
-- pyyaml
-- pytest (for development)
-
-## License
-
-[License details to be added] 
+## References
+- Chapman-Richards growth model literature
+- Jimenez, J.A. & Lugo, A.E. (1985, 1987)
+- Zanvo et al. (2023) for allometric equations
+- Field data from Nigeria mangrove sites (2024)

configs/{params.yaml → declining_increment.yaml}

@@ -1,76 +1,80 @@
+growth_model: declining_increment
+continuous_growth: true
+
 species:
   - name: "species_A"  # Rhizophora spp.
     planting_density: 4000  # trees/ha
     mortality_rates:
-      year_1: 15  # %
-      year_2: 12  # %
-      year_3: 6   # %
+      year_1: 20  # %
+      year_2: 15  # %
+      year_3: 7   # %
       year_4: 5   # %
       year_5: 5   # %
-      subsequent: 1.5  # %
-    chapman_richards:
+      year_6: 4   # %
+      year_7: 3   # %
+      year_8: 2   # %
+      year_9: 1.5 # %
+      year_10: 1  # %
+      subsequent: 1  # % (years 11-30)
+    declining_increment:
       dbh:
-        a: 11.07  # asymptotic maximum DBH (cm)
-        b: 0.25   # growth rate parameter (yr^-1)
-        c: 2.1    # shape parameter
+        r0: 0.58  # cm/yr
+        T_m: 40  # years
       height:
-        a: 12.0   # asymptotic maximum height (m)
-        b: 0.25   # growth rate parameter (yr^-1)
-        c: 2.1    # shape parameter
+        r0: 0.43  # m/yr
+        T_m: 40  # years
     allometry:
-      equation: "Zanvo et al. 2023: Total = 1.938 × (DBH² H)^0.67628"
+      equation: "Zanvo et al. 2023: Total = 2.0738 × (DBH² H)^0.67628"
     initial_values:
-      dbh: 1.0    # cm
-      height: 0.5  # m
+      dbh: 0.1    # cm
+      height: 0.2  # m
 
   - name: "species_B"  # Avicennia germinans
     planting_density: 444  # trees/ha
     mortality_rates:
-      year_1: 15  # %
-      year_2: 12  # %
-      year_3: 6   # %
+      year_1: 20  # %
+      year_2: 15  # %
+      year_3: 7   # %
       year_4: 5   # %
       year_5: 5   # %
-      subsequent: 1.5  # %
-    chapman_richards:
+      year_6: 4   # %
+      year_7: 3   # %
+      year_8: 2   # %
+      year_9: 1.5 # %
+      year_10: 1  # %
+      subsequent: 1  # % (years 11-30)
+    declining_increment:
       dbh:
-        a: 17.0   # asymptotic maximum DBH (cm)
-        b: 0.25   # growth rate parameter (yr^-1)
-        c: 2.1    # shape parameter
+        r0: 0.91  # cm/yr
+        T_m: 40
       height:
-        a: 8.8    # asymptotic maximum height (m)
-        b: 0.25   # growth rate parameter (yr^-1)
-        c: 2.1    # shape parameter
+        r0: 0.46  # m/yr
+        T_m: 40
     allometry:
-      equation: "Zanvo et al. 2023: Total = 1.486 × (DBH² H)^0.55864"
+      equation: "Zanvo et al. 2023: Total = 1.5595 × (DBH² H)^0.55864"
     initial_values:
-      dbh: 1.0    # cm
-      height: 0.5  # m
+      dbh: 0.1    # cm
+      height: 0.2  # m
 
 project:
   duration_years: 30  # matches project_years
   planting_schedule:
-    year_1: 2500  # ha
-    year_2: 2000  # ha
-    year_3: 0  # ha
-    year_4: 0  # ha
-    year_5: 0  # ha
+    year_1: 500
+    year_2: 1500
+    year_3: 1500
+    year_4: 1500
+    year_5: 0
 
 carbon:
   biomass_to_carbon: 0.47  # IPCC default
   carbon_to_co2: 3.67     # molecular weight ratio
-  buffer_percentage: 15    # risk buffer
+  buffer_percentage: 20    # risk buffer
   leakage_percentage: 0    # no leakage assumed
   baseline_emissions: 0    # tCO2/ha/year
-
-# Note: b parameter is calculated using b = ln(c) / Tm
-# For both species:
-# - Tm = 6 years
-# - c = 2.1
-# Therefore: b = ln(2.1) / 6 ≈ 0.35 
+  soil_carbon_per_ha_per_year: 1.0  # tCO2/ha/year, added annually
 
 scenarios:
   area: 1000.0  # ha
-  dbh_range: [1.0, 25.0]  # cm
-  height_range: [0.5, 15.0]  # m
+  dbh_range: [0.1, 25.0]  # cm
+  height_range: [0.2, 15.0]  # m
   growth_rate_factor: 1.0 

dashboard/app.py

@@ -1,404 +1,693 @@
 """
-Gradio dashboard for exploring ER model results and parameters.
+Mangrove ER Model Dashboard
+
+- Each tab corresponds to a different model/config.
+- To add a new model, add its config to configs/ and add an entry to MODEL_CONFIGS below.
+- Each tab shows the results for the base config (no parameter editing).
 """
-from pathlib import Path
-from typing import Dict, List, Tuple
 
 import gradio as gr
-import matplotlib.pyplot as plt
-import numpy as np
+from pathlib import Path
+import sys
+import os
+
+# Add the parent directory to Python path so we can import er_model_core
+sys.path.insert(0, str(Path(__file__).parent.parent))
+
+from er_model_core.er_model import ERModel
 import pandas as pd
+import numpy as np
+import plotly.graph_objs as go
+from plotly.subplots import make_subplots
+import warnings
+import traceback
+from er_model_core.allometry import calculate_biomass
+from dashboard.plots import check_complex, create_growth_increment_plots, create_all_plots
+from dashboard.tables import create_survival_table, create_biomass_table, create_project_results_table
 import yaml
+import tempfile
+from dashboard.inputs import get_model_inputs
 
-from src.er_model import ERModel
+MODEL_CONFIGS = {
+    "Declining Increment": "configs/declining_increment.yaml"
+}
 
+# Add a mapping for display names
+SPECIES_DISPLAY_NAMES = {
+    'species_A': 'Rhizophora spp.',
+    'species_B': 'Avicennia germinans'
+}
+
+def to_native(obj):
+    if isinstance(obj, pd.DataFrame):
+        obj = obj.applymap(lambda x: to_native(x))
+        obj.index = obj.index.map(to_native)
+        obj.columns = obj.columns.map(to_native)
+        return obj
+    elif isinstance(obj, pd.Series):
+        obj = obj.map(to_native)
+        obj.index = obj.index.map(to_native)
+        return obj
+    elif isinstance(obj, dict):
+        return {to_native(k): to_native(v) for k, v in obj.items()}
+    elif isinstance(obj, list):
+        return [to_native(x) for x in obj]
+    elif isinstance(obj, tuple):
+        return tuple(to_native(x) for x in obj)
+    elif isinstance(obj, (np.generic, np.ndarray)):
+        return obj.item() if obj.size == 1 else obj.tolist()
+    else:
+        return obj
 
-def load_config(config_path: Path = Path("configs/params.yaml")) -> dict:
-    """Load and parse YAML configuration."""
+# Helper to update planting schedule in a config file
+def update_planting_schedule(config_path, year_areas):
     with open(config_path) as f:
-        return yaml.safe_load(f)
-
-
-def run_model(config: dict) -> Tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame]:
-    """Run model with current parameters."""
-    # Save temporary config
-    temp_config = Path("configs/temp_config.yaml")
-    with open(temp_config, "w") as f:
-        yaml.dump(config, f)
-    
-    # Run model and cleanup
-    model = ERModel(temp_config)
-    results, species_results, scenario_results = model.run()
-    temp_config.unlink()
-    
-    return results, species_results, scenario_results
-
-
-def format_number(x):
-    """Format number with comma separators, handling NaN and floating points."""
-    if pd.isna(x):
-        return "0"
-    if isinstance(x, (int, float)):
-        if isinstance(x, int) or float(x).is_integer():
-            return format(int(x), ',')
-        return format(round(float(x), 2), ',')
-    return str(x)
-
-
-def create_plots_and_tables(
-    results: pd.DataFrame,
-    species_results: pd.DataFrame,
-    scenario_results: pd.DataFrame
-) -> Tuple[plt.Figure, plt.Figure, plt.Figure, str, pd.DataFrame, pd.DataFrame]:
-    """Create carbon sequestration plots, summary, and tables."""
-    # Clear any existing plots
-    plt.close('all')
-    
-    # Time series plot
-    fig1 = plt.figure(figsize=(10, 6))
-    ax1 = fig1.add_subplot(111)
-    results.plot(
-        x="year",
-        y=["gross_carbon", "net_carbon"],
-        ax=ax1,
-        title="Carbon Sequestration Over Time"
-    )
-    ax1.set_xlabel("Year")
-    ax1.set_ylabel("Carbon (tCO2)")
-    ax1.legend(["Gross Carbon", "Net Carbon"])
-    # Format y-axis with comma separator
-    ax1.yaxis.set_major_formatter(plt.FuncFormatter(lambda x, p: format_number(x)))
-    
-    # Milestone years plot
-    milestones = [5, 10, 15, 20, 25, 30]
-    milestone_data = results[results["year"].isin(milestones)]
-    
-    fig2 = plt.figure(figsize=(10, 6))
-    ax2 = fig2.add_subplot(111)
-    width = 0.35
-    x = np.arange(len(milestones))
-    ax2.bar(x - width/2, milestone_data["gross_carbon"], width, label="Gross Carbon")
-    ax2.bar(x + width/2, milestone_data["net_carbon"], width, label="Net Carbon")
-    ax2.set_xticks(x)
-    ax2.set_xticklabels(milestones)
-    ax2.set_xlabel("Year")
-    ax2.set_ylabel("Carbon (tCO2)")
-    ax2.set_title("Carbon Sequestration at Milestone Years")
-    ax2.legend()
-    # Format y-axis with comma separator
-    ax2.yaxis.set_major_formatter(plt.FuncFormatter(lambda x, p: format_number(x)))
-
-    # New biomass per tree plot
-    fig3 = plt.figure(figsize=(10, 6))
-    ax3 = fig3.add_subplot(111)
-    
-    # Load config to get planting densities
-    config = load_config()
-    rhiz_density = config["species"][0]["planting_density"]  # Rhizophora density
-    avic_density = config["species"][1]["planting_density"]  # Avicennia density
-    
-    # Extract biomass per tree data from species_results
-    years = species_results.index
-    
-    # Get species column names (first two columns after 'Year' should be the species)
-    species_cols = [col for col in species_results.columns if 'tCO2' in col]
-    if len(species_cols) >= 2:
-        rhiz_col = species_cols[0]  # First species column
-        avic_col = species_cols[1]  # Second species column
-        
-        # Calculate biomass per tree using actual planting densities
-        rhiz_biomass = species_results[rhiz_col].astype(float) / (results["area_year_1"].iloc[0] * rhiz_density)
-        avic_biomass = species_results[avic_col].astype(float) / (results["area_year_1"].iloc[0] * avic_density)
-        
-        ax3.plot(years, rhiz_biomass, 'o-', color='#1f77b4', label='Rhizophora spp. - Total Biomass')
-        ax3.plot(years, avic_biomass, 'o-', color='#ff7f0e', label='Avicennia germinans - Total Biomass')
-        
-        ax3.set_xlabel("Year since planting")
-        ax3.set_ylabel("Biomass (tonnes per tree)")
-        ax3.set_title("Total Biomass per Tree")
-        ax3.grid(True)
-        ax3.legend()
-
-    # Calculate total area from planting schedule
-    total_area = sum(float(v) for k, v in results.iloc[-1].items() if k.startswith("area_year"))
-    
-    # Get final year carbon values
+        config = yaml.safe_load(f)
+    for i, area in enumerate(year_areas, 1):
+        config["project"]["planting_schedule"][f"year_{i}"] = area
+    return config
+
+def create_survival_table(model):
+    """
+    Returns a DataFrame with years as rows and columns for each species and total surviving trees.
+    Uses model.species_metrics for all values.
+    """
+    df = model.species_metrics.copy()
+    # Pivot to wide format: years as rows, species as columns
+    surv = df.pivot(index="Year", columns="Species", values="Surviving Trees")
+    # Rename columns to display names
+    surv.columns = [SPECIES_DISPLAY_NAMES.get(sp, sp) for sp in surv.columns]
+    # Add total surviving trees column
+    surv["Total Surviving Trees"] = surv.sum(axis=1)
+    # Format numbers
+    for col in surv.columns:
+        surv[col] = surv[col].apply(lambda x: f"{x:,.0f}")
+    surv = surv.reset_index()
+    return surv
+
+def run_model_from_config(config):
+    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+        yaml.dump(config, tmp)
+        tmp_path = tmp.name
+    model = ERModel(Path(tmp_path))
+    results, species_results = model.run()
+    plots = create_all_plots(results, species_results, config)
+    summary = create_summary(results, species_results, config, model)
+    survival_table = create_survival_table(model)
+    # Use modularized table formatting
+    results_fmt = create_project_results_table(results)
+    species_results_fmt = species_results.copy()
+    for col in species_results_fmt.columns:
+        if species_results_fmt[col].dtype in [float, int]:
+            species_results_fmt[col] = species_results_fmt[col].apply(lambda x: f"{x:,.2f}" if isinstance(x, float) else f"{x:,}")
+    return (*plots, summary, results_fmt, species_results_fmt, survival_table)
+
+def run_model(config_path):
+    with open(config_path) as f:
+        config = yaml.safe_load(f)
+    model = ERModel(Path(config_path))
+    results, species_results = model.run()
+    plots = create_all_plots(results, species_results, config)
+    summary = create_summary(results, species_results, config, model)
+    survival_table = create_survival_table(model)
+    # Use modularized table formatting
+    results_fmt = create_project_results_table(results)
+    species_results_fmt = species_results.copy()
+    for col in species_results_fmt.columns:
+        if species_results_fmt[col].dtype in [float, int]:
+            species_results_fmt[col] = species_results_fmt[col].apply(lambda x: f"{x:,.2f}" if isinstance(x, float) else f"{x:,}")
+    return (*plots, summary, results_fmt, species_results_fmt, survival_table)
+
+def create_summary(results, species_results, config, model):
+    total_area = sum(config["project"]["planting_schedule"].values())
     final_gross = results["gross_carbon"].iloc[-1]
-    final_net = results["net_carbon"].iloc[-1]
-    
-    # Calculate per hectare metrics
-    if total_area > 0:
-        gross_carbon_per_ha = final_gross / total_area
-        net_carbon_per_ha = final_net / total_area
+    final_buffer = results["buffer_carbon"].iloc[-1]
+    # Sum soil carbon over all years
+    total_soil_carbon = results["soil_carbon"].sum()
+    final_gross_soil = final_gross + total_soil_carbon
+    final_buffer_soil = final_buffer + (total_soil_carbon * (1 - config["carbon"]["buffer_percentage"] / 100))
+    years = len(results)
+    gross_carbon_per_ha = final_gross / total_area if total_area > 0 else 0
+    buffer_carbon_per_ha = final_buffer / total_area if total_area > 0 else 0
+    gross_carbon_per_ha_soil = final_gross_soil / total_area if total_area > 0 else 0
+    buffer_carbon_per_ha_soil = final_buffer_soil / total_area if total_area > 0 else 0
+    annual_buffer_per_ha = buffer_carbon_per_ha / years if years > 0 else 0
+    annual_buffer_per_ha_soil = buffer_carbon_per_ha_soil / years if years > 0 else 0
+    soil_carbon_val = config["carbon"].get("soil_carbon_per_ha_per_year", 0)
+    # Milestone years
+    def get_val(col, idx):
+        return results[col].iloc[idx] if len(results) > idx else 0
+    year_5_gross = get_val("gross_carbon", 4)
+    year_5_buffer = get_val("buffer_carbon", 4)
+    year_10_gross = get_val("gross_carbon", 9)
+    year_10_buffer = get_val("buffer_carbon", 9)
+    year_20_gross = get_val("gross_carbon", 19)
+    year_20_buffer = get_val("buffer_carbon", 19)
+    avg_annual_gross = final_gross / years if years > 0 else 0
+    avg_annual_buffer = final_buffer / years if years > 0 else 0
+    avg_annual_gross_soil = final_gross_soil / years if years > 0 else 0
+    avg_annual_buffer_soil = final_buffer_soil / years if years > 0 else 0
+    # Surviving trees at milestones
+    milestones = [5, 10, 20, 30]
+    survival_lines = []
+    for m in milestones:
+        if m <= years:
+            surv = model.calculate_total_surviving_trees(m)
+            surv_str = ", ".join([f"{SPECIES_DISPLAY_NAMES.get(k, k)}: {v:,.0f}" for k, v in surv.items()])
+            survival_lines.append(f"Year {m} surviving trees: {surv_str}")
+    # Per hectare (final year)
+    surv_30 = model.calculate_total_surviving_trees(years)
+    per_ha_lines = []
+    for k, v in surv_30.items():
+        per_ha = v / total_area if total_area > 0 else 0
+        per_ha_lines.append(f"{SPECIES_DISPLAY_NAMES.get(k, k)}: {per_ha:,.0f} trees/ha")
+    summary = (
+        f"Project Overview:\n----------------\nDuration: {years} years\nTotal Area Planted: {total_area:,.0f} ha\nBuffer Pool: {config['carbon']['buffer_percentage']}%\n\n"
+        f"Soil Carbon:\n-----------\nSoil carbon added: {soil_carbon_val} tCO2/ha/year\nTotal soil carbon added over project: {total_soil_carbon:,.0f} tCO2\n\n"
+        f"Carbon Sequestration (Biomass Only):\n-------------------\nTotal Gross Carbon (Year {years}): {final_gross:,.0f} tCO2\nTotal Buffer Carbon (Year {years}): {final_buffer:,.0f} tCO2\nAverage Annual Gross: {avg_annual_gross:,.0f} tCO2/yr\nAverage Annual Buffer: {avg_annual_buffer:,.0f} tCO2/yr\n\n"
+        f"Carbon Sequestration (With Soil Carbon):\n-------------------\nTotal Gross Carbon (Year {years}): {final_gross_soil:,.0f} tCO2\nTotal Buffer Carbon (Year {years}): {final_buffer_soil:,.0f} tCO2\nAverage Annual Gross: {avg_annual_gross_soil:,.0f} tCO2/yr\nAverage Annual Buffer: {avg_annual_buffer_soil:,.0f} tCO2/yr\n\n"
+        f"Milestone Years (Biomass Only):\n--------------\nYear 5 Gross: {year_5_gross:,.0f} tCO2\nYear 5 Buffer: {year_5_buffer:,.0f} tCO2\nYear 10 Gross: {year_10_gross:,.0f} tCO2\nYear 10 Buffer: {year_10_buffer:,.0f} tCO2\nYear 20 Gross: {year_20_gross:,.0f} tCO2\nYear 20 Buffer: {year_20_buffer:,.0f} tCO2\n\n"
+        "Surviving Trees (Milestones):\n----------------------------\n"
+        + "\n".join(survival_lines)
+        + "\n\nPer Hectare Surviving Trees (Final Year):\n----------------------------------------\n"
+        + "\n".join(per_ha_lines)
+        + f"\n\nPer Hectare Metrics (Year {years}):\n-----------------------------\nGross Carbon per ha (biomass only): {gross_carbon_per_ha:,.0f} tCO2/ha\nBuffer Carbon per ha (biomass only): {buffer_carbon_per_ha:,.0f} tCO2/ha\nGross Carbon per ha (with soil): {gross_carbon_per_ha_soil:,.0f} tCO2/ha\nBuffer Carbon per ha (with soil): {buffer_carbon_per_ha_soil:,.0f} tCO2/ha\nAverage Annual Buffer per ha (biomass only): {annual_buffer_per_ha:,.2f} tCO2/ha/yr\nAverage Annual Buffer per ha (with soil): {annual_buffer_per_ha_soil:,.2f} tCO2/ha/yr\n"
+    )
+    return summary
+
+def hex_to_rgba(hex_color, alpha):
+    hex_color = hex_color.lstrip('#')
+    r, g, b = tuple(int(hex_color[i:i+2], 16) for i in (0, 2, 4))
+    return f'rgba({r},{g},{b},{alpha})'
+
+def get_all_model_inputs(config):
+    # Use get_model_inputs to get all widgets and defaults
+    model_inputs, model_defaults = get_model_inputs(config)
+    # Unpack widgets for layout
+    planting_density_A = model_inputs['planting_density_A']
+    r0_dbh_A = model_inputs['r0_dbh_A']
+    tm_dbh_A = model_inputs['tm_dbh_A']
+    r0_height_A = model_inputs['r0_height_A']
+    tm_height_A = model_inputs['tm_height_A']
+    planting_density_B = model_inputs['planting_density_B']
+    r0_dbh_B = model_inputs['r0_dbh_B']
+    tm_dbh_B = model_inputs['tm_dbh_B']
+    r0_height_B = model_inputs['r0_height_B']
+    tm_height_B = model_inputs['tm_height_B']
+    mort_vars = [model_inputs[f'year_{yr}'] for yr in range(1, 11)]
+    mort_sub = model_inputs['subsequent']
+    project_duration = model_inputs['project_duration']
+    buffer_pct = model_inputs['buffer_pct']
+    soil_carbon = model_inputs['soil_carbon']
+    return (planting_density_A, r0_dbh_A, tm_dbh_A, r0_height_A, tm_height_A,
+            planting_density_B, r0_dbh_B, tm_dbh_B, r0_height_B, tm_height_B,
+            mort_vars, mort_sub, project_duration, buffer_pct, soil_carbon, model_defaults)
+
+# Gradio UI code starts here
+with gr.Blocks() as demo:
+    gr.Markdown("""
+    <div style='display: flex; align-items: center; gap: 16px;'>
+        <img src='dashboard/sequestr_logo.png' alt='Sequestr Logo' style='height:48px;'>
+        <div>
+            <b>This tool is provided by <a href='https://www.sequestr.io/' target='_blank'>Sequestr</a>.</b><br>
+            <span style='color: #d9534f;'>Please do not share this page.</span>
+        </div>
+    </div>
+    """, elem_id="sequestr-header")
+    gr.Markdown(r"""
+# Mangrove Emissions Reduction (ER) Model & Dashboard
+
+Welcome to the Mangrove ER Model Dashboard, a tool designed by [Sequestr](https://www.sequestr.io/) for simulating and analyzing carbon sequestration in our mangrove restoration project. This dashboard allows users to explore the "Declining Increment" growth model, adjust key project parameters, and visualize the resulting carbon and emission reduction outcomes over time.
+""")
+    with gr.Accordion("How the Model Works", open=False):
+        gr.Markdown(r"""
+## How the Model Works
+
+The core of this dashboard is an emissions reduction (ER) model that projects mangrove growth, biomass accumulation, and carbon sequestration based on user-defined parameters and a chosen growth model. Currently, the active model is the **Declining Increment** model.
+
+### Declining Increment Growth Model
+
+This model assumes that the annual growth increment (for both Diameter at Breast Height - DBH, and Height) declines linearly over time, eventually reaching zero. The total size at any given age is the sum of all non-negative annual increments up to that age. This ensures that tree size never decreases.
+
+**1. Discrete Declining Increment:** The size at a given age $t$ (in years) is calculated by summing annual increments:
+
+$$
+\mathrm{size}(t) = \mathrm{initial\_size} + \sum_{i=1}^{t} \mathrm{increment}(i)
+$$
+
+Where the increment for year $i$ is:
+
+$$
+\mathrm{increment}(i) = \max \left(0, r_0 \cdot \left(1 - \frac{i - 0.5}{T_m}\right) \right)
+$$
+
+- $\mathrm{initial\_size}$: The initial DBH (cm) or Height (m) at planting (age 0).
+- $r_0$: The initial annual growth increment (cm/year for DBH, m/year for Height).
+- $T_m$: The time (in years) at which the annual growth increment becomes zero.
+
+**2. Continuous Declining Increment (Analytical Formula):** For a smoother representation, the continuous version uses an analytical formula:
+
+$$
+\mathrm{size}(t) = \mathrm{initial\_size} + r_0 \cdot \left(t - \frac{t^2}{2T_m}\right), \quad 0 \leq t \leq T_m
+$$
+
+If $t > T_m$, then $\mathrm{size}(t) = \mathrm{size}(T_m)$. This version may yield slightly different results from the discrete sum, especially for short time periods or small $T_m$ values. The dashboard configuration specifies which version is active.
+
+### Survival and Mortality
+The number of surviving trees in a cohort is calculated annually. Mortality can be specified using either:
+
+1. **Annual Mortality Rates:** A percentage mortality rate applied each year, potentially varying for the first few years and then stabilizing.
+2. **DBH-Dependent Mortality (Optional, if configured):** Mortality rate $m$ is a function of tree DBH:
+
+$$
+m = m_\mathrm{ref} \cdot \left( \frac{\mathrm{DBH}_\mathrm{ref}}{\mathrm{DBH}} \right)^p
+$$
+
+- $m_\mathrm{ref}$: Reference mortality rate at $\mathrm{DBH}_\mathrm{ref}$.
+- $\mathrm{DBH}_\mathrm{ref}$: Reference DBH.
+- $p$: Exponent controlling the sensitivity of mortality to DBH. The calculated $m$ is capped (e.g., between 0 and 0.99).
+""")
+    with gr.Accordion("How Emission Reductions (ERs) are Calculated", open=False):
+        gr.Markdown(r"""
+## How Emission Reductions (ERs) are Calculated
+
+The model follows these steps to estimate net CO2 emission reductions:
+
+**1. Biomass Calculation:** Total above-ground and below-ground biomass per tree is calculated using allometric equations, which typically relate DBH and Height to biomass. For example, using equations from Zanvo et al. (2023):
+- For *Rhizophora racemosa* (Species A):
+  $$
+  \mathrm{Biomass}_{\mathrm{total}} = 2.0738 \times (\mathrm{DBH}^2 \cdot H)^{0.67628}
+  $$
+- For *Avicennia germinans* (Species B):
+  $$
+  \mathrm{Biomass}_{\mathrm{total}} = 1.5595 \times (\mathrm{DBH}^2 \cdot H)^{0.55864}
+  $$
+  *(Note: The specific equations are defined in the model configuration.)*
+
+**2. Carbon Stock Calculation:** The total carbon stock in living biomass per hectare is calculated:
+
+$$
+\mathrm{Carbon\ Stock\ (tC/ha)} = \frac{\mathrm{Biomass\ per\ tree\ (kg)} \times \mathrm{Surviving\ Trees\ per\ ha} \times \mathrm{Biomass\ to\ Carbon\ ratio}}{1000}
+$$
+
+- The Biomass to Carbon ratio (e.g., 0.47) converts biomass to carbon mass.
+
+**3. Gross CO2 Sequestration:** The carbon stock is then converted to tons of CO2 equivalent:
+
+$$
+\mathrm{Gross\ CO2eq\ (tCO2/ha)} = \mathrm{Carbon\ Stock\ (tC/ha)} \times \mathrm{Carbon\ to\ CO2\ ratio}
+$$
+
+- The Carbon to CO2 ratio (typically 3.67) is based on molecular weights.
+
+**4. Soil Carbon (Optional):** If configured, annual soil carbon sequestration is added:
+
+$$
+\mathrm{Soil\ Carbon\ CO2eq\ (tCO2/ha/yr)} = \mathrm{User\ Defined\ Value\ (e.g.\ 1.0\ tCO2/ha/yr)}
+$$
+
+This is added to the gross CO2 sequestration from biomass.
+
+**5. Net Emission Reductions (ERs):** Adjustments are made to the gross CO2 sequestration (including soil carbon, if applicable) to determine net ERs eligible for crediting:
+
+$$
+\mathrm{Net\ ERs} = (\mathrm{Gross\ CO2eq}_{\mathrm{biomass + soil}}) \times (1 - \mathrm{Buffer\ \%}) - (\mathrm{Gross\ CO2eq}_{\mathrm{biomass + soil}} \times \mathrm{Leakage\ \%}) - (\mathrm{Baseline\ Emissions\ per\ ha} \times \mathrm{Area})
+$$
+
+- **Buffer Pool**: A percentage deduction to account for risks like project failure or natural disturbances.
+- **Leakage**: Emissions occurring outside the project boundary due to project activities (often assumed to be 0% for mangrove projects if activities are self-contained).
+- **Baseline Emissions**: Emissions that would have occurred in the absence of the project (e.g., from degrading land).
+
+This dashboard visualizes these values annually over the project duration, providing insights into the project's carbon sequestration potential.
+""")
+    # --- Key Input Params Card ---
+    config = yaml.safe_load(open(MODEL_CONFIGS["Declining Increment"]))
+    species_md = "| Species | Planting Density (trees/ha) | r0 (DBH, cm/yr) | Tm (DBH, yrs) | r0 (Height, m/yr) | Tm (Height, yrs) |\n|---|---|---|---|---|---|\n"
+    for sp in config["species"]:
+        species_md += f"| {sp['name']} | {sp['planting_density']} | {sp['declining_increment']['dbh']['r0']} | {sp['declining_increment']['dbh']['T_m']} | {sp['declining_increment']['height']['r0']} | {sp['declining_increment']['height']['T_m']} |\n"
+    # Project details with proper line breaks
+    project_md = (
+        f"**Project Duration:** {config['project']['duration_years']} years  \n"
+        f"**Buffer %:** {config['carbon']['buffer_percentage']}  \n"
+        f"**Soil Carbon per ha per year:** {config['carbon']['soil_carbon_per_ha_per_year']} tCO2"
+    )
+    # Mortality table
+    mortality_schedule = config['project'].get('annual_mortality_schedule', {})
+    mortality_md = "| Year | Mortality Rate (%) |\n|---|---|\n"
+    if mortality_schedule:
+        for year, rate in mortality_schedule.items():
+            mortality_md += f"| {year} | {rate} |\n"
     else:
-        gross_carbon_per_ha = 0
-        net_carbon_per_ha = 0
-    
-    # Create summary text with formatted numbers
-    summary = f"""
-Results Summary:
----------------
-Total Years: {len(results)}
-Total Area Planted: {format_number(total_area)} ha
-Total Gross Carbon: {format_number(final_gross)} tCO2
-Total Net Carbon: {format_number(final_net)} tCO2
-
-Per Hectare Metrics (Year 30):
------------------------------
-Gross Carbon per ha: {format_number(gross_carbon_per_ha)} tCO2/ha
-Net Carbon per ha: {format_number(net_carbon_per_ha)} tCO2/ha
-"""
-    
-    # Format species results table with comma separators
-    species_table = species_results.copy()
-    species_table.index.name = "Year"
-    for col in species_table.select_dtypes(include=[np.number]).columns:
-        species_table[col] = species_table[col].apply(format_number)
-    
-    # Format scenario comparison table with comma separators
-    scenario_table = scenario_results.copy()
-    scenario_table.index.name = "Scenario"
-    for col in scenario_table.select_dtypes(include=[np.number]).columns:
-        scenario_table[col] = scenario_table[col].apply(format_number)
-    
-    return fig1, fig2, fig3, summary, species_table, scenario_table
-
-
-def update_model(
-    rhiz_density: float,
-    rhiz_mort_1: float,
-    rhiz_mort_2: float,
-    rhiz_mort_3: float,
-    rhiz_mort_4: float,
-    rhiz_mort_5: float,
-    avic_density: float,
-    avic_mort_1: float,
-    avic_mort_2: float,
-    avic_mort_3: float,
-    avic_mort_4: float,
-    avic_mort_5: float,
-    year_1_area: float,
-    year_2_area: float,
-    year_3_area: float,
-    year_4_area: float,
-    year_5_area: float,
-    buffer_percentage: float,
-    scenario_area: float = 1000.0,
-    dbh_range_min: float = 1.0,
-    dbh_range_max: float = 25.0,
-    height_range_min: float = 0.5,
-    height_range_max: float = 15.0,
-    growth_rate_factor: float = 1.0,
-) -> Tuple[plt.Figure, plt.Figure, plt.Figure, str, pd.DataFrame, pd.DataFrame]:
-    """Update model with new parameters and return plots and tables."""
-    # Load base config
-    config = load_config()
-    
-    # Update Rhizophora parameters
-    config["species"][0]["planting_density"] = rhiz_density
-    config["species"][0]["mortality_rates"].update({
-        "year_1": rhiz_mort_1,
-        "year_2": rhiz_mort_2,
-        "year_3": rhiz_mort_3,
-        "year_4": rhiz_mort_4,
-        "year_5": rhiz_mort_5
-    })
-    
-    # Update Avicennia parameters
-    config["species"][1]["planting_density"] = avic_density
-    config["species"][1]["mortality_rates"].update({
-        "year_1": avic_mort_1,
-        "year_2": avic_mort_2,
-        "year_3": avic_mort_3,
-        "year_4": avic_mort_4,
-        "year_5": avic_mort_5
-    })
-    
-    # Update planting schedule
-    config["project"]["planting_schedule"].update({
-        "year_1": year_1_area,
-        "year_2": year_2_area,
-        "year_3": year_3_area,
-        "year_4": year_4_area,
-        "year_5": year_5_area
-    })
-    
-    # Update carbon parameters
-    config["carbon"]["buffer_percentage"] = buffer_percentage
-    
-    # Add scenario parameters
-    config["scenarios"] = {
-        "area": scenario_area,
-        "dbh_range": [dbh_range_min, dbh_range_max],
-        "height_range": [height_range_min, height_range_max],
-        "growth_rate_factor": growth_rate_factor
-    }
-    
-    # Run model and create plots/tables
-    results, species_results, scenario_results = run_model(config)
-    return create_plots_and_tables(results, species_results, scenario_results)
-
-
-def main():
-    """Launch the Gradio interface."""
-    config = load_config()
-    
-    # Extract default values
-    rhiz_defaults = config["species"][0]
-    avic_defaults = config["species"][1]
-    planting_defaults = config["project"]["planting_schedule"]
-    carbon_defaults = config["carbon"]
-    scenario_defaults = config.get("scenarios", {
-        "area": 1000.0,
-        "dbh_range": [1.0, 25.0],
-        "height_range": [0.5, 15.0],
-        "growth_rate_factor": 1.0
-    })
-    
-    # Create interface
-    with gr.Blocks(title="ER Model Dashboard", theme=gr.themes.Soft()) as interface:
-        gr.Markdown("# ER Model Dashboard")
-        gr.Markdown("Explore carbon sequestration results and test different parameters")
-        
-        with gr.Row():
-            with gr.Column():
-                # Rhizophora parameters
-                gr.Markdown("### Rhizophora Parameters")
-                rhiz_density = gr.Slider(
-                    minimum=100,
-                    maximum=10000,
-                    value=rhiz_defaults["planting_density"],
-                    label="Planting Density (trees/ha)"
-                )
-                rhiz_mort_1 = gr.Slider(0, 100, rhiz_defaults["mortality_rates"]["year_1"], label="Year 1 Mortality (%)")
-                rhiz_mort_2 = gr.Slider(0, 100, rhiz_defaults["mortality_rates"]["year_2"], label="Year 2 Mortality (%)")
-                rhiz_mort_3 = gr.Slider(0, 100, rhiz_defaults["mortality_rates"]["year_3"], label="Year 3 Mortality (%)")
-                rhiz_mort_4 = gr.Slider(0, 100, rhiz_defaults["mortality_rates"]["year_4"], label="Year 4 Mortality (%)")
-                rhiz_mort_5 = gr.Slider(0, 100, rhiz_defaults["mortality_rates"]["year_5"], label="Year 5 Mortality (%)")
+        mortality_md += "| (no data) | (no data) |\n"
+    # Carbon details
+    carbon = config["carbon"]
+    carbon_md = (
+        f"**Biomass to Carbon:** {carbon['biomass_to_carbon']}  \n"
+        f"**Carbon to CO2:** {carbon['carbon_to_co2']}  \n"
+        f"**Buffer %:** {carbon['buffer_percentage']}  \n"
+        f"**Leakage %:** {carbon['leakage_percentage']}  \n"
+        f"**Baseline Emissions:** {carbon['baseline_emissions']} tCO2/ha/year  \n"
+        f"**Soil Carbon per ha per year:** {carbon['soil_carbon_per_ha_per_year']} tCO2"
+    )
+    gr.Markdown("""
+### Key Model Inputs
+""" + species_md)
+    gr.Markdown("""
+### Carbon Details
+""" + carbon_md)
+    with gr.Tabs():
+        with gr.Tab("Declining Increment"):
+            config = yaml.safe_load(open(MODEL_CONFIGS["Declining Increment"]))
+            use_continuous = config.get('continuous_growth', False)
+            if use_continuous:
+                gr.Markdown("**Mode:** Continuous declining increment growth")
+            else:
+                gr.Markdown("**Mode:** Discrete declining increment growth (sum of annual increments)")
+            gr.Markdown("## Model Inputs")
+            with gr.Row():
+                # Column 1: Planting Schedule (as gr.Number)
+                with gr.Column(scale=1, min_width=220):
+                    with gr.Accordion("Planting Schedule (ha per year)", open=True):
+                        year_1_default = config['project']['planting_schedule']['year_1']
+                        year_2_default = config['project']['planting_schedule']['year_2']
+                        year_3_default = config['project']['planting_schedule']['year_3']
+                        year_4_default = config['project']['planting_schedule']['year_4']
+                        year_5_default = config['project']['planting_schedule']['year_5']
+                        year_1 = gr.Number(value=year_1_default, label="Area to be planted in year 1 (ha)", interactive=True)
+                        year_2 = gr.Number(value=year_2_default, label="Area to be planted in year 2 (ha)", interactive=True)
+                        year_3 = gr.Number(value=year_3_default, label="Area to be planted in year 3 (ha)", interactive=True)
+                        year_4 = gr.Number(value=year_4_default, label="Area to be planted in year 4 (ha)", interactive=True)
+                        year_5 = gr.Number(value=year_5_default, label="Area to be planted in year 5 (ha)", interactive=True)
+                # Column 2: Planting Density (as gr.Slider)
+                with gr.Column(scale=1, min_width=220):
+                    with gr.Accordion("Planting Density", open=True):
+                        pdA = config['species'][0]['planting_density']
+                        pdB = config['species'][1]['planting_density']
+                        planting_density_A = gr.Slider(minimum=pdA*0.75, maximum=pdA*1.25, value=pdA, step=1, label="Initial planting density for Species A")
+                        planting_density_B = gr.Slider(minimum=pdB*0.75, maximum=pdB*1.25, value=pdB, step=1, label="Initial planting density for Species B")
+                # Columns 3-4: Mortality Rates (as gr.Slider, split across two half-width columns)
+                with gr.Column(scale=2, min_width=300):
+                    with gr.Accordion("Mortality Rates (%)", open=True):
+                        mort_defaults = [config['species'][0]['mortality_rates'][f'year_{yr}'] for yr in range(1, 11)]
+                        mort_sub_default = config['species'][0]['mortality_rates']['subsequent']
+                        with gr.Row():
+                            with gr.Column():
+                                mort_vars_1 = [gr.Slider(minimum=val*0.75, maximum=val*1.25, value=val, step=0.1, label=f"Annual mortality rate for year {i+1}") for i, val in enumerate(mort_defaults[:6])]
+                            with gr.Column():
+                                mort_vars_2 = [gr.Slider(minimum=val*0.75, maximum=val*1.25, value=val, step=0.1, label=f"Annual mortality rate for year {i+7}") for i, val in enumerate(mort_defaults[6:])]
+                                mort_sub = gr.Slider(minimum=mort_sub_default*0.75, maximum=mort_sub_default*1.25, value=mort_sub_default, step=0.1, label="Annual mortality rate for years 11+")
+                        mort_vars = mort_vars_1 + mort_vars_2
+            # After the mortality columns, add:
+            soil_carbon_default = config['carbon']['soil_carbon_per_ha_per_year']
+
+            soil_carbon = gr.Number(value=soil_carbon_default, label="Soil Carbon per ha per year (tCO2)", interactive=True, info="Annual soil carbon addition per hectare")
+
+            # Add a Reset to Defaults button
+            reset_btn = gr.Button("Reset to Defaults", variant="secondary")
+
+            def reset_to_defaults():
+                # Reload config and return all default values in widget order
+                config = yaml.safe_load(open(MODEL_CONFIGS["Declining Increment"]))
+                pdA = config['species'][0]['planting_density']
+                pdB = config['species'][1]['planting_density']
+                mort_defaults = [config['species'][0]['mortality_rates'][f'year_{yr}'] for yr in range(1, 11)]
+                mort_sub_default = config['species'][0]['mortality_rates']['subsequent']
+                return [
+                    config['project']['planting_schedule']['year_1'],
+                    config['project']['planting_schedule']['year_2'],
+                    config['project']['planting_schedule']['year_3'],
+                    config['project']['planting_schedule']['year_4'],
+                    config['project']['planting_schedule']['year_5'],
+                    pdA, pdB,
+                    *mort_defaults, mort_sub_default,
+                    config['carbon']['soil_carbon_per_ha_per_year']
+                ]
+            reset_btn.click(
+                reset_to_defaults,
+                inputs=[],
+                outputs=[year_1, year_2, year_3, year_4, year_5,
+                         planting_density_A, planting_density_B,
+                         *mort_vars, mort_sub, soil_carbon]
+            )
+
+            update_btn = gr.Button("Update Results", variant="primary")
+            with gr.Row():
+                carbon_plot = gr.Plot()
+                biomass_plot = gr.Plot()
+                annual_plot = gr.Plot()
+            growth_plot = gr.Plot(label="Growth & Increment Curves (DBH/Height)")
+            summary_box = gr.Textbox(label="Summary", lines=12)
+            # --- Tabbed tables section ---
+            with gr.Tabs():
+                with gr.Tab("Project Results (Annual)"):
+                    results_box = gr.Dataframe(label="Project Results (Annual)")
+                    download_results_file = gr.File(label="Download CSV", visible=False)
+                    download_results_btn = gr.Button("📥 Download CSV", variant="secondary")
+                with gr.Tab("Species Results (Annual)"):
+                    species_box = gr.Dataframe(label="Species Results (Annual)")
+                    download_species_file = gr.File(label="Download CSV", visible=False)
+                    download_species_btn = gr.Button("📥 Download CSV", variant="secondary")
+                with gr.Tab("Surviving Trees Table"):
+                    survival_box = gr.Dataframe(label="Surviving Trees Table")
+                    download_survival_file = gr.File(label="Download CSV", visible=False)
+                    download_survival_btn = gr.Button("📥 Download CSV", variant="secondary")
+                with gr.Tab("Biomass Table"):
+                    biomass_debug_table = gr.Dataframe(label="Biomass Table (inputs & outputs per year)")
+                    download_biomass_file = gr.File(label="Download CSV", visible=False)
+                    download_biomass_btn = gr.Button("📥 Download CSV", variant="secondary")
+            # --- End tabbed tables ---
+            # Update the update_declining_increment callback to use these new inputs
+            def update_declining_increment(y1, y2, y3, y4, y5,
+                                           planting_density_A, planting_density_B,
+                                           *mortality_inputs_and_soil_carbon):
+                # Unpack the last one from the end
+                *mortality_inputs, soil_carbon = mortality_inputs_and_soil_carbon
+                config = update_planting_schedule(MODEL_CONFIGS["Declining Increment"], [y1, y2, y3, y4, y5])
+                # Update species_A params
+                config['species'][0]['planting_density'] = planting_density_A
+                # Update species_B params
+                config['species'][1]['planting_density'] = planting_density_B
+                # Update mortality for BOTH species from the same UI values (years 1-10 + subsequent)
+                for i in [0, 1]:
+                    for yr_idx, yr in enumerate(range(1, 11)):
+                        config['species'][i]['mortality_rates'][f'year_{yr}'] = mortality_inputs[yr_idx]
+                    config['species'][i]['mortality_rates']['subsequent'] = mortality_inputs[10]
+                # Update project/carbon params
+                config['carbon']['soil_carbon_per_ha_per_year'] = soil_carbon
+                with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+                    yaml.dump(config, tmp)
+                    tmp_path = tmp.name
+                model = ERModel(Path(tmp_path))
+                results, species_results = model.run()
+                plots = create_all_plots(results, species_results, config)
+                summary = create_summary(results, species_results, config, model)
+                survival_table = create_survival_table(model)
+                growth_fig = create_growth_increment_plots(config, model_type="declining_increment")
+                biomass_debug_df = create_biomass_table(config)
+                results_fmt = create_project_results_table(results)
+                species_results_fmt = species_results.copy()
+                for col in species_results_fmt.columns:
+                    if species_results_fmt[col].dtype in [float, int]:
+                        species_results_fmt[col] = species_results_fmt[col].apply(lambda x: f"{x:,.2f}" if isinstance(x, float) else f"{x:,}")
+                # --- Ensure all DataFrames are native types ---
+                results_fmt = to_native(results_fmt)
+                species_results_fmt = to_native(species_results_fmt)
+                survival_table = to_native(survival_table)
+                biomass_debug_df = to_native(biomass_debug_df)
+                
+                # --- Ensure all other values are native types ---
+                if hasattr(summary, 'item'):
+                    summary = summary.item()
+                return (plots[0], plots[1], plots[2], growth_fig, summary, results_fmt, species_results_fmt, survival_table, biomass_debug_df)
+            update_btn.click(
+                update_declining_increment,
+                inputs=[year_1, year_2, year_3, year_4, year_5,
+                        planting_density_A, planting_density_B,
+                        *mort_vars, mort_sub, soil_carbon],
+                outputs=[carbon_plot, biomass_plot, annual_plot, growth_plot, summary_box, results_box, species_box, survival_box, biomass_debug_table]
+            )
+
+            # Download functionality for tables
+            def create_download_file(df, filename):
+                """Create a temporary CSV file for download"""
+                import tempfile
+                import os
+                
+                # Create a temporary file
+                temp_dir = tempfile.gettempdir()
+                filepath = os.path.join(temp_dir, filename)
+                
+                # Convert formatted strings back to numbers where possible
+                if df is not None and not df.empty:
+                    df_clean = df.copy()
+                    for col in df_clean.columns:
+                        if df_clean[col].dtype == 'object' and col.lower() != 'year':
+                            try:
+                                # Try to convert formatted strings back to numbers
+                                df_clean[col] = df_clean[col].astype(str).str.replace(',', '').astype(float)
+                            except:
+                                pass  # Keep as string if conversion fails
+                    
+                    df_clean.to_csv(filepath, index=False)
+                    return filepath
+                return None
+            
+            def download_results_csv(y1, y2, y3, y4, y5, planting_density_A, planting_density_B, *mortality_inputs_and_soil_carbon):
+                """Generate and download project results CSV"""
+                try:
+                    # Recreate the results using current parameters
+                    *mortality_inputs, soil_carbon = mortality_inputs_and_soil_carbon
+                    config = update_planting_schedule(MODEL_CONFIGS["Declining Increment"], [y1, y2, y3, y4, y5])
+                    config['species'][0]['planting_density'] = planting_density_A
+                    config['species'][1]['planting_density'] = planting_density_B
+                    for i in [0, 1]:
+                        for yr_idx, yr in enumerate(range(1, 11)):
+                            config['species'][i]['mortality_rates'][f'year_{yr}'] = mortality_inputs[yr_idx]
+                        config['species'][i]['mortality_rates']['subsequent'] = mortality_inputs[10]
+                    config['carbon']['soil_carbon_per_ha_per_year'] = soil_carbon
+                    
+                    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+                        yaml.dump(config, tmp)
+                        tmp_path = tmp.name
+                    model = ERModel(Path(tmp_path))
+                    results, _ = model.run()
+                    
+                    filepath = create_download_file(results, "project_results.csv")
+                    return gr.File(value=filepath, visible=True) if filepath else gr.File(visible=False)
+                except Exception as e:
+                    print(f"Error creating download: {e}")
+                    return gr.File(visible=False)
                 
-            with gr.Column():
-                # Avicennia parameters
-                gr.Markdown("### Avicennia Parameters")
-                avic_density = gr.Slider(
-                    minimum=100,
-                    maximum=10000,
-                    value=avic_defaults["planting_density"],
-                    label="Planting Density (trees/ha)"
-                )
-                avic_mort_1 = gr.Slider(0, 100, avic_defaults["mortality_rates"]["year_1"], label="Year 1 Mortality (%)")
-                avic_mort_2 = gr.Slider(0, 100, avic_defaults["mortality_rates"]["year_2"], label="Year 2 Mortality (%)")
-                avic_mort_3 = gr.Slider(0, 100, avic_defaults["mortality_rates"]["year_3"], label="Year 3 Mortality (%)")
-                avic_mort_4 = gr.Slider(0, 100, avic_defaults["mortality_rates"]["year_4"], label="Year 4 Mortality (%)")
-                avic_mort_5 = gr.Slider(0, 100, avic_defaults["mortality_rates"]["year_5"], label="Year 5 Mortality (%)")
-        
-        with gr.Row():
-            with gr.Column():
-                # Planting schedule
-                gr.Markdown("### Planting Schedule")
-                year_1_area = gr.Slider(0, 1000, planting_defaults["year_1"], label="Year 1 Area (ha)")
-                year_2_area = gr.Slider(0, 1000, planting_defaults["year_2"], label="Year 2 Area (ha)")
-                year_3_area = gr.Slider(0, 1000, planting_defaults["year_3"], label="Year 3 Area (ha)")
-                year_4_area = gr.Slider(0, 1000, planting_defaults["year_4"], label="Year 4 Area (ha)")
-                year_5_area = gr.Slider(0, 1000, planting_defaults["year_5"], label="Year 5 Area (ha)")
+            def download_species_csv(y1, y2, y3, y4, y5, planting_density_A, planting_density_B, *mortality_inputs_and_soil_carbon):
+                """Generate and download species results CSV"""
+                try:
+                    *mortality_inputs, soil_carbon = mortality_inputs_and_soil_carbon
+                    config = update_planting_schedule(MODEL_CONFIGS["Declining Increment"], [y1, y2, y3, y4, y5])
+                    config['species'][0]['planting_density'] = planting_density_A
+                    config['species'][1]['planting_density'] = planting_density_B
+                    for i in [0, 1]:
+                        for yr_idx, yr in enumerate(range(1, 11)):
+                            config['species'][i]['mortality_rates'][f'year_{yr}'] = mortality_inputs[yr_idx]
+                        config['species'][i]['mortality_rates']['subsequent'] = mortality_inputs[10]
+                    config['carbon']['soil_carbon_per_ha_per_year'] = soil_carbon
+                    
+                    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+                        yaml.dump(config, tmp)
+                        tmp_path = tmp.name
+                    model = ERModel(Path(tmp_path))
+                    _, species_results = model.run()
+                    
+                    filepath = create_download_file(species_results, "species_results.csv")
+                    return gr.File(value=filepath, visible=True) if filepath else gr.File(visible=False)
+                except Exception as e:
+                    print(f"Error creating download: {e}")
+                    return gr.File(visible=False)
                 
-            with gr.Column():
-                # Carbon parameters
-                gr.Markdown("### Carbon Parameters")
-                buffer_percentage = gr.Slider(
-                    minimum=0,
-                    maximum=30,
-                    value=carbon_defaults["buffer_percentage"],
-                    label="Buffer Percentage (%)"
-                )
+            def download_survival_csv(y1, y2, y3, y4, y5, planting_density_A, planting_density_B, *mortality_inputs_and_soil_carbon):
+                """Generate and download survival table CSV"""
+                try:
+                    *mortality_inputs, soil_carbon = mortality_inputs_and_soil_carbon
+                    config = update_planting_schedule(MODEL_CONFIGS["Declining Increment"], [y1, y2, y3, y4, y5])
+                    config['species'][0]['planting_density'] = planting_density_A
+                    config['species'][1]['planting_density'] = planting_density_B
+                    for i in [0, 1]:
+                        for yr_idx, yr in enumerate(range(1, 11)):
+                            config['species'][i]['mortality_rates'][f'year_{yr}'] = mortality_inputs[yr_idx]
+                        config['species'][i]['mortality_rates']['subsequent'] = mortality_inputs[10]
+                    config['carbon']['soil_carbon_per_ha_per_year'] = soil_carbon
+                    
+                    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+                        yaml.dump(config, tmp)
+                        tmp_path = tmp.name
+                    model = ERModel(Path(tmp_path))
+                    model.run()
+                    survival_table = create_survival_table(model)
+                    
+                    # Convert survival table to raw numbers
+                    df_raw = model.species_metrics.copy()
+                    surv_raw = df_raw.pivot(index="Year", columns="Species", values="Surviving Trees")
+                    surv_raw.columns = [SPECIES_DISPLAY_NAMES.get(sp, sp) for sp in surv_raw.columns]
+                    surv_raw["Total Surviving Trees"] = surv_raw.sum(axis=1)
+                    surv_raw = surv_raw.reset_index()
+                    
+                    filepath = create_download_file(surv_raw, "surviving_trees.csv")
+                    return gr.File(value=filepath, visible=True) if filepath else gr.File(visible=False)
+                except Exception as e:
+                    print(f"Error creating download: {e}")
+                    return gr.File(visible=False)
                 
-                # Scenario parameters
-                gr.Markdown("### Scenario Parameters")
-                scenario_area = gr.Slider(
-                    minimum=100,
-                    maximum=2000,
-                    value=scenario_defaults["area"],
-                    label="Scenario Area (ha)"
-                )
-                dbh_range_min = gr.Slider(
-                    minimum=0.5,
-                    maximum=10,
-                    value=scenario_defaults["dbh_range"][0],
-                    label="Min DBH (cm)"
-                )
-                dbh_range_max = gr.Slider(
-                    minimum=10,
-                    maximum=30,
-                    value=scenario_defaults["dbh_range"][1],
-                    label="Max DBH (cm)"
-                )
-                height_range_min = gr.Slider(
-                    minimum=0.5,
-                    maximum=5,
-                    value=scenario_defaults["height_range"][0],
-                    label="Min Height (m)"
-                )
-                height_range_max = gr.Slider(
-                    minimum=5,
-                    maximum=20,
-                    value=scenario_defaults["height_range"][1],
-                    label="Max Height (m)"
-                )
-                growth_rate_factor = gr.Slider(
-                    minimum=0.5,
-                    maximum=2.0,
-                    value=scenario_defaults["growth_rate_factor"],
-                    label="Growth Rate Factor"
-                )
-        
-        # Outputs
-        with gr.Row():
-            plot1 = gr.Plot(label="Carbon Sequestration Over Time")
-            plot2 = gr.Plot(label="Milestone Years")
-        
-        with gr.Row():
-            plot3 = gr.Plot(label="Total Biomass per Tree")
-        
-        summary = gr.Textbox(label="Summary", lines=10)
-        
-        with gr.Row():
-            species_table = gr.Dataframe(
-                label="Per Species Carbon Change & Cumulative Hectares",
-                headers=["Year", "Rhizophora tCO2", "Avicennia tCO2", "Total tCO2", "Cumulative ha", "tCO2/ha"]
+            def download_biomass_csv(y1, y2, y3, y4, y5, planting_density_A, planting_density_B, *mortality_inputs_and_soil_carbon):
+                """Generate and download biomass table CSV"""
+                try:
+                    *mortality_inputs, soil_carbon = mortality_inputs_and_soil_carbon
+                    config = update_planting_schedule(MODEL_CONFIGS["Declining Increment"], [y1, y2, y3, y4, y5])
+                    config['species'][0]['planting_density'] = planting_density_A
+                    config['species'][1]['planting_density'] = planting_density_B
+                    for i in [0, 1]:
+                        for yr_idx, yr in enumerate(range(1, 11)):
+                            config['species'][i]['mortality_rates'][f'year_{yr}'] = mortality_inputs[yr_idx]
+                        config['species'][i]['mortality_rates']['subsequent'] = mortality_inputs[10]
+                    config['carbon']['soil_carbon_per_ha_per_year'] = soil_carbon
+                    
+                    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+                        yaml.dump(config, tmp)
+                        tmp_path = tmp.name
+                    model = ERModel(Path(tmp_path))
+                    model.run()
+                    
+                    filepath = create_download_file(model.species_metrics, "biomass_table.csv")
+                    return gr.File(value=filepath, visible=True) if filepath else gr.File(visible=False)
+                except Exception as e:
+                    print(f"Error creating download: {e}")
+                    return gr.File(visible=False)
+            
+            # Connect download buttons to functions
+            download_results_btn.click(
+                download_results_csv,
+                inputs=[year_1, year_2, year_3, year_4, year_5,
+                        planting_density_A, planting_density_B,
+                        *mort_vars, mort_sub, soil_carbon],
+                outputs=[download_results_file]
             )
-            scenario_table = gr.Dataframe(
-                label="Scenario Comparison",
-                headers=["Scenario", "DBH Range", "Height Range", "Growth Rate", "Year 5 tCO2", "Year 10 tCO2", "Year 15 tCO2", "Year 30 tCO2"]
+            
+            download_species_btn.click(
+                download_species_csv,
+                inputs=[year_1, year_2, year_3, year_4, year_5,
+                        planting_density_A, planting_density_B,
+                        *mort_vars, mort_sub, soil_carbon],
+                outputs=[download_species_file]
             )
-        
-        # Update function
-        inputs = [
-            rhiz_density, rhiz_mort_1, rhiz_mort_2, rhiz_mort_3, rhiz_mort_4, rhiz_mort_5,
-            avic_density, avic_mort_1, avic_mort_2, avic_mort_3, avic_mort_4, avic_mort_5,
-            year_1_area, year_2_area, year_3_area, year_4_area, year_5_area,
-            buffer_percentage, scenario_area, dbh_range_min, dbh_range_max,
-            height_range_min, height_range_max, growth_rate_factor
-        ]
-        
-        outputs = [plot1, plot2, plot3, summary, species_table, scenario_table]
-        
-        # Initialize with default values
-        interface.load(
-            fn=update_model,
-            inputs=inputs,
-            outputs=outputs,
-        )
-        
-        # Update on any change
-        for inp in inputs:
-            inp.change(
-                fn=update_model,
-                inputs=inputs,
-                outputs=outputs
+            
+            download_survival_btn.click(
+                download_survival_csv,
+                inputs=[year_1, year_2, year_3, year_4, year_5,
+                        planting_density_A, planting_density_B,
+                        *mort_vars, mort_sub, soil_carbon],
+                outputs=[download_survival_file]
             )
-    
-    # Launch interface
-    interface.launch(share=False)
-
-
-if __name__ == "__main__":
-    main() 
+            
+            download_biomass_btn.click(
+                download_biomass_csv,
+                inputs=[year_1, year_2, year_3, year_4, year_5,
+                        planting_density_A, planting_density_B,
+                        *mort_vars, mort_sub, soil_carbon],
+                outputs=[download_biomass_file]
+            )
+            
+            # Show initial results
+            c, b, a, g, summary, r, s, surv, biomass_debug_df = update_declining_increment(
+                year_1_default, year_2_default, year_3_default, year_4_default, year_5_default,
+                pdA, pdB,
+                *mort_defaults, mort_sub_default,
+                soil_carbon_default
+            )
+            carbon_plot.value = c
+            biomass_plot.value = b
+            annual_plot.value = a
+            growth_plot.value = g
+            summary_box.value = summary
+            results_box.value = r
+            species_box.value = s
+            survival_box.value = surv
+            biomass_debug_table.value = biomass_debug_df
+    demo.launch(share=True) 

dashboard/inputs.py

@@ -0,0 +1,56 @@
+import gradio as gr
+
+def get_model_inputs(config):
+    spA = config["species"][0]
+    spB = config["species"][1]
+    # Species A
+    planting_density_A = gr.Number(value=spA['planting_density'], label=f"{spA['name']} Planting Density (trees/ha)")
+    r0_dbh_A = gr.Number(value=spA['declining_increment']['dbh']['r0'], label=f"{spA['name']} r0 (DBH, cm/yr)")
+    tm_dbh_A = gr.Number(value=spA['declining_increment']['dbh']['T_m'], label=f"{spA['name']} Tm (DBH, yrs)")
+    r0_height_A = gr.Number(value=spA['declining_increment']['height']['r0'], label=f"{spA['name']} r0 (Height, m/yr)")
+    tm_height_A = gr.Number(value=spA['declining_increment']['height']['T_m'], label=f"{spA['name']} Tm (Height, yrs)")
+    # Species B
+    planting_density_B = gr.Number(value=spB['planting_density'], label=f"{spB['name']} Planting Density (trees/ha)")
+    r0_dbh_B = gr.Number(value=spB['declining_increment']['dbh']['r0'], label=f"{spB['name']} r0 (DBH, cm/yr)")
+    tm_dbh_B = gr.Number(value=spB['declining_increment']['dbh']['T_m'], label=f"{spB['name']} Tm (DBH, yrs)")
+    r0_height_B = gr.Number(value=spB['declining_increment']['height']['r0'], label=f"{spB['name']} r0 (Height, m/yr)")
+    tm_height_B = gr.Number(value=spB['declining_increment']['height']['T_m'], label=f"{spB['name']} Tm (Height, yrs)")
+    # Mortality (one set, mirrored, years 1-10 + subsequent)
+    mort_inputs = {}
+    mort_defaults = []
+    for yr in range(1, 11):
+        key = f"year_{yr}"
+        mort_inputs[key] = gr.Number(value=spA['mortality_rates'][key], label=f"Year {yr} Mortality (%)")
+        mort_defaults.append(spA['mortality_rates'][key])
+    mort_inputs['subsequent'] = gr.Number(value=spA['mortality_rates']['subsequent'], label="Subsequent Years Mortality (%)")
+    mort_defaults.append(spA['mortality_rates']['subsequent'])
+    # Project/Carbon
+    project_duration = gr.Number(value=config['project']['duration_years'], label="Project Duration (years)")
+    buffer_pct = gr.Number(value=config['carbon']['buffer_percentage'], label="Buffer %")
+    soil_carbon = gr.Number(value=config['carbon']['soil_carbon_per_ha_per_year'], label="Soil Carbon per ha per year (tCO2)")
+    # Return as dict and list of default values (for initial call)
+    components = {
+        'planting_density_A': planting_density_A,
+        'r0_dbh_A': r0_dbh_A,
+        'tm_dbh_A': tm_dbh_A,
+        'r0_height_A': r0_height_A,
+        'tm_height_A': tm_height_A,
+        'planting_density_B': planting_density_B,
+        'r0_dbh_B': r0_dbh_B,
+        'tm_dbh_B': tm_dbh_B,
+        'r0_height_B': r0_height_B,
+        'tm_height_B': tm_height_B,
+        **mort_inputs,
+        'project_duration': project_duration,
+        'buffer_pct': buffer_pct,
+        'soil_carbon': soil_carbon
+    }
+    defaults = [
+        spA['planting_density'], spA['declining_increment']['dbh']['r0'], spA['declining_increment']['dbh']['T_m'],
+        spA['declining_increment']['height']['r0'], spA['declining_increment']['height']['T_m'],
+        spB['planting_density'], spB['declining_increment']['dbh']['r0'], spB['declining_increment']['dbh']['T_m'],
+        spB['declining_increment']['height']['r0'], spB['declining_increment']['height']['T_m'],
+        *mort_defaults,
+        config['project']['duration_years'], config['carbon']['buffer_percentage'], config['carbon']['soil_carbon_per_ha_per_year']
+    ]
+    return components, defaults 

dashboard/plots/__init__.py

@@ -0,0 +1,134 @@
+"""
+Plotting functions for the Mangrove ER Dashboard.
+"""
+# Imports will be updated as needed when moving functions from app.py 
+import plotly.graph_objects as go
+from plotly.subplots import make_subplots
+import numpy as np
+import warnings
+import traceback
+from er_model_core.allometry import calculate_biomass
+from er_model_core.er_model import ERModel
+from pathlib import Path
+from er_model_core.growth_models.chapman_richards import chapman_richards_growth
+from er_model_core.growth_models.linear import linear_growth, linear_plateau_growth
+from er_model_core.growth_models.declining_increment import declining_increment_growth, continuous_declining_increment_growth
+
+# Helper for checking complex numbers (copied from app.py)
+def check_complex(arr, label):
+    if np.iscomplexobj(arr):
+        complex_indices = np.where(np.iscomplex(arr))[0]
+        warnings.warn(f"[ERROR] Complex values in {label}: indices={complex_indices}, values={arr[complex_indices]}")
+        traceback.print_stack()
+
+def create_growth_increment_plots(config, model_type=None):
+    """
+    Create a 2x2 grid of growth and increment plots for DBH and Height using Plotly.
+    model_type: 'chapman_richards', 'linear', 'linear_plateau', or 'declining_increment'
+    """
+    if model_type is None:
+        model_type = config.get('growth_model', 'Unknown Model')
+    N = config["project"]["duration_years"]
+    ages = np.arange(0, N + 1)  # 0 to N
+    ages_inc = np.arange(1, N + 1)  # 1 to N
+    fig = make_subplots(rows=2, cols=2, subplot_titles=("DBH Growth", "HEIGHT Growth", "DBH Annual Increment", "HEIGHT Annual Increment"))
+    from er_model_core.er_model import ERModel
+    SPECIES_DISPLAY_NAMES = {
+        'species_A': 'Rhizophora spp.',
+        'species_B': 'Avicennia germinans'
+    }
+    use_continuous = config.get('continuous_growth', False)
+    for sp in config["species"]:
+        name = SPECIES_DISPLAY_NAMES.get(sp["name"], sp["name"])
+        initial_dbh = sp["initial_values"]["dbh"]
+        initial_height = sp["initial_values"]["height"]
+        if model_type == "linear":
+            dbh = [linear_growth(t, sp["linear"]["dbh"], initial_dbh) for t in ages]
+            height = [linear_growth(t, sp["linear"]["height"], initial_height) for t in ages]
+        elif model_type == "linear_plateau":
+            dbh = [linear_plateau_growth(t, sp["linear_plateau"]["dbh"], initial_dbh) for t in ages]
+            height = [linear_plateau_growth(t, sp["linear_plateau"]["height"], initial_height) for t in ages]
+        elif model_type == "declining_increment":
+            if use_continuous:
+                dbh = [continuous_declining_increment_growth(t, sp["declining_increment"]["dbh"], initial_dbh) for t in ages]
+                height = [continuous_declining_increment_growth(t, sp["declining_increment"]["height"], initial_height) for t in ages]
+            else:
+                dbh = [declining_increment_growth(t, sp["declining_increment"]["dbh"], initial_dbh) for t in ages]
+                height = [declining_increment_growth(t, sp["declining_increment"]["height"], initial_height) for t in ages]
+        else:
+            dbh = [chapman_richards_growth(t, sp["chapman_richards"]["dbh"], initial_dbh) for t in ages]
+            height = [chapman_richards_growth(t, sp["chapman_richards"]["height"], initial_height) for t in ages]
+        dbh = np.array(dbh)
+        height = np.array(height)
+        check_complex(dbh, f"{name} DBH (growth)")
+        check_complex(height, f"{name} Height (growth)")
+        dbh_inc = dbh[1:] - dbh[:-1]
+        height_inc = height[1:] - height[:-1]
+        check_complex(dbh_inc, f"{name} DBH Δ (increment)")
+        check_complex(height_inc, f"{name} Height Δ (increment)")
+        fig.add_trace(go.Scatter(x=ages, y=dbh, mode='lines', name=f"{name} DBH", legendgroup=name, line=dict(width=2)), row=1, col=1)
+        fig.add_trace(go.Scatter(x=ages, y=height, mode='lines', name=f"{name} Height", legendgroup=name, line=dict(width=2)), row=1, col=2)
+        fig.add_trace(go.Scatter(x=ages_inc, y=dbh_inc, mode='lines', name=f"{name} DBH Δ", legendgroup=name, showlegend=False, line=dict(width=2)), row=2, col=1)
+        fig.add_trace(go.Scatter(x=ages_inc, y=height_inc, mode='lines', name=f"{name} Height Δ", legendgroup=name, showlegend=False, line=dict(width=2)), row=2, col=2)
+    fig.update_layout(height=700, width=900, title_text="Growth and Increment Curves", hovermode="x unified")
+    fig.update_xaxes(title_text="Age (years)", row=1, col=1)
+    fig.update_xaxes(title_text="Age (years)", row=1, col=2)
+    fig.update_xaxes(title_text="Age (years)", row=2, col=1)
+    fig.update_xaxes(title_text="Age (years)", row=2, col=2)
+    fig.update_yaxes(title_text="Size (cm)", row=1, col=1)
+    fig.update_yaxes(title_text="Size (m)", row=1, col=2)
+    fig.update_yaxes(title_text="Annual increment (cm/year)", row=2, col=1)
+    fig.update_yaxes(title_text="Annual increment (m/year)", row=2, col=2)
+    return fig
+
+
+def create_all_plots(results, species_results, config):
+    # Diagnostic: Check for complex numbers in results DataFrame columns used for plotting
+    for col in ["gross_carbon", "buffer_carbon"]:
+        arr = np.array(results[col])
+        check_complex(arr, f"results['{col}']")
+    # 1. Carbon curve (gross and buffer, with and without soil)
+    fig1 = go.Figure()
+    fig1.add_trace(go.Scatter(x=results["year"], y=results["gross_carbon"], mode="lines+markers", name="Gross Carbon", line=dict(width=2, color="blue")))
+    fig1.add_trace(go.Scatter(x=results["year"], y=results["buffer_carbon"], mode="lines+markers", name="Buffer Carbon", line=dict(width=2, color="red")))
+    if "gross_carbon_with_soil" in results.columns:
+        # For backward compatibility, but we now sum soil carbon in summary, so recalc here for plot
+        soil_cumsum = results["soil_carbon"].cumsum()
+        gross_with_soil = results["gross_carbon"] + soil_cumsum
+        buffer_with_soil = results["buffer_carbon"] + (soil_cumsum * (1 - config["carbon"]["buffer_percentage"] / 100))
+        fig1.add_trace(go.Scatter(x=results["year"], y=gross_with_soil, mode="lines+markers", name="Gross Carbon (with Soil)", line=dict(width=2, dash="dot", color="green")))
+        fig1.add_trace(go.Scatter(x=results["year"], y=buffer_with_soil, mode="lines+markers", name="Buffer Carbon (with Soil)", line=dict(width=2, dash="dot", color="#43c6ac")))
+    fig1.update_layout(title="Carbon Sequestration Over Time", xaxis_title="Year", yaxis_title="Carbon (tCO2)", hovermode="x unified", template="plotly_white")
+
+    # 2. Annual buffer carbon (ERs)
+    annual_ers = results["buffer_carbon"].diff().fillna(results["buffer_carbon"].iloc[0])
+    arr = np.array(annual_ers)
+    check_complex(arr, "annual_ers")
+    fig2 = go.Figure()
+    fig2.add_trace(go.Bar(x=results["year"], y=annual_ers, name="Annual Buffer ERs", marker_color="#2ecc71", opacity=0.7))
+    fig2.update_layout(title="Annual Emission Reductions (Post Buffer)", xaxis_title="Year", yaxis_title="Annual Buffer ERs (tCO2)", hovermode="x unified", template="plotly_white")
+
+    # 3. Biomass per tree for all species (refactored to use model.species_metrics)
+    years = results["year"]
+    fig3 = go.Figure()
+    import tempfile, yaml
+    from er_model_core.er_model import ERModel
+    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+        yaml.dump(config, tmp)
+        tmp_path = tmp.name
+    model = ERModel(Path(tmp_path))
+    model.run()
+    df = model.species_metrics.copy()
+    SPECIES_DISPLAY_NAMES = {
+        'species_A': 'Rhizophora spp.',
+        'species_B': 'Avicennia germinans'
+    }
+    for sp in config["species"]:
+        name = SPECIES_DISPLAY_NAMES.get(sp["name"], sp["name"])
+        # Get per-tree biomass for this species
+        sp_biomass = df[df["Species"] == sp["name"]].set_index("Year")["Biomass per Tree (kg)"]
+        # Ensure correct order and fill missing years with NaN
+        sp_biomass = sp_biomass.reindex(years)
+        fig3.add_trace(go.Scatter(x=years, y=sp_biomass, mode="lines+markers", name=name, line=dict(width=2)))
+    fig3.update_layout(title="Total Biomass per Tree", xaxis_title="Year since planting", yaxis_title="Biomass (kg per tree)", hovermode="x unified", template="plotly_white")
+    return (fig1, fig2, fig3) 

dashboard/tables/__init__.py

@@ -0,0 +1,79 @@
+"""
+Table formatting functions for the Mangrove ER Dashboard.
+"""
+# Imports will be updated as needed when moving functions from app.py 
+import pandas as pd
+import numpy as np
+from er_model_core.er_model import ERModel
+from pathlib import Path
+import yaml
+
+# Helper for pretty table titles
+def pretty_table_title(title):
+    return f"<span style='font-size:1.3em; font-weight:bold'>{title}</span>"
+
+# Table: Biomass Table (was debug table)
+def create_biomass_table(config):
+    """
+    Returns a DataFrame with years as rows and columns for each species:
+    - Surviving Trees
+    - DBH (cm)
+    - Height (m)
+    - Biomass per Tree (kg)
+    - Total Biomass (kg)
+    Uses model.species_metrics for all values.
+    """
+    import tempfile
+    with tempfile.NamedTemporaryFile(mode="w", suffix=".yaml", delete=False) as tmp:
+        yaml.dump(config, tmp)
+        tmp_path = tmp.name
+    model = ERModel(Path(tmp_path))
+    model.run()
+    df = model.species_metrics.copy()
+    # Pivot to multi-index columns: (Species, Metric)
+    df = df.pivot(index="Year", columns="Species")
+    # Flatten columns
+    SPECIES_DISPLAY_NAMES = {
+        'species_A': 'Rhizophora spp.',
+        'species_B': 'Avicennia germinans'
+    }
+    df.columns = [f"{SPECIES_DISPLAY_NAMES.get(sp, sp)}, {metric}" for metric, sp in df.columns]
+    # Format numbers
+    for col in df.columns:
+        if "Surviving Trees" in col:
+            df[col] = df[col].apply(lambda x: f"{x:,.0f}")
+        elif "Biomass" in col or "DBH" in col or "Height" in col:
+            df[col] = df[col].apply(lambda x: f"{x:,.2f}")
+    return df
+
+# Table: Surviving Trees Table
+def create_survival_table(model):
+    """
+    Returns a DataFrame with years as rows and columns for each species and total surviving trees.
+    Uses model.species_metrics for all values.
+    """
+    SPECIES_DISPLAY_NAMES = {
+        'species_A': 'Rhizophora spp.',
+        'species_B': 'Avicennia germinans'
+    }
+    df = model.species_metrics.copy()
+    # Pivot to wide format: years as rows, species as columns
+    surv = df.pivot(index="Year", columns="Species", values="Surviving Trees")
+    # Rename columns to display names
+    surv.columns = [SPECIES_DISPLAY_NAMES.get(sp, sp) for sp in surv.columns]
+    # Add total surviving trees column
+    surv["Total Surviving Trees"] = surv.sum(axis=1)
+    # Format numbers
+    for col in surv.columns:
+        surv[col] = surv[col].apply(lambda x: f"{x:,.0f}")
+    surv = surv.reset_index()
+    return surv
+
+def create_project_results_table(results: pd.DataFrame) -> pd.DataFrame:
+    df = results.copy()
+    # No need to cumsum soil_carbon; model now provides cumulative values
+    # Optionally format numbers for display
+    for col in df.columns:
+        if df[col].dtype in [float, int]:
+            df[col] = df[col].apply(lambda x: f"{x:,.2f}" if isinstance(x, float) else f"{x:,}")
+    return df 

environment.yml

@@ -1,20 +1,16 @@
 name: er-model
 channels:
   - conda-forge
+  - defaults
 dependencies:
   - python=3.10
-  - numpy=1.24.3
-  - pandas=2.0.3
-  - scipy=1.12.0
-  - matplotlib=3.8.2
-  - pyyaml=6.0.1
-  - click=8.1.7
-  - pip=24.0
-  - pytest=7.4.4
-  - black=24.1.1
-  - isort=5.13.2
-  - pre-commit=3.6.0
+  - pip
   - pip:
-    - gradio==3.50.2
-    - fastapi==0.95.2
-    - pydantic==1.10.13 
+    - click>=8.1.7
+    - gradio>=4.0.0
+    - matplotlib>=3.8.2
+    - numpy>=1.26.3
+    - pandas>=2.2.0
+    - plotly>=5.0.0
+    - pyyaml>=6.0.1
+    - scipy>=1.12.0 

{src → er_model_core}/__init__.py

@@ -2,9 +2,9 @@
 ER Model package for calculating carbon sequestration in mangrove projects.
 """
 
-from .allometry import calculate_biomass
-from .er_model import ERModel, Species, ProjectConfig, CarbonConfig
-from .metrics import calculate_carbon
+from er_model_core.allometry import calculate_biomass
+from er_model_core.er_model import ERModel, Species, ProjectConfig, CarbonConfig
+from er_model_core.metrics import calculate_carbon
 
 __all__ = [
     "ERModel",

{src → er_model_core}/allometry.py

@@ -17,15 +17,21 @@ def calculate_biomass(dbh: float, height: float, species_name: str, params: Dict
     Returns:
         Total tree biomass (above + below ground) in kg
     """
+    base = dbh**2 * height
+    if base <= 0:
+        print(f"[ERROR] Non-positive base in calculate_biomass: dbh={dbh}, height={height}, base={base}, species={species_name}")
+
     # Use Zanvo et al. 2023 equations that include both DBH and height
-    if species_name == "species_A":  # Rhizophora
-        # Total = 1.938 × (DBH² H)^0.67628
-        total_biomass = 1.938 * (dbh**2 * height)**0.67628
-    elif species_name == "species_B":  # Avicennia
-        # Total = 1.486 × (DBH² H)^0.55864
-        total_biomass = 1.486 * (dbh**2 * height)**0.55864
+    # Use Komiyama et al. 2005 generic mangrove equation for the else statement, the generic for all other mangroves
+    # Note that Komiyama generic equation is used with the average density value of 0.71 g/cm3 found across various mangrove species
+    if species_name == "species_A":  # Rhizophora racemosa
+        # Total = 2.0738 × (DBH² H)^0.67628
+        total_biomass = 2.0738 * (dbh**2 * height)**0.67628
+    elif species_name == "species_B":  # Avicennia germinans
+        # Total = 1.5595 × (DBH² H)^0.55864
+        total_biomass = 1.5595 * (dbh**2 * height)**0.55864
     else:
         # Use a conservative generic equation if species not recognized
-        total_biomass = 1.712 * (dbh**2 * height)**0.61746
+        total_biomass = 0.251 * 0.71 * (dbh**2.46)
     
     return total_biomass 

er_model_core/config_loader.py

@@ -0,0 +1,46 @@
+"""
+Configuration loading and parsing for the ER Model.
+"""
+from dataclasses import dataclass
+from pathlib import Path
+from typing import Dict, List, Tuple, Optional
+import yaml
+
+# Assuming these dataclasses are defined in er_model.er_model or a shared types module.
+# For now, let's duplicate them here or import them if they are easily accessible.
+# Ideally, they'd be in a place like 'er_model.types'
+# For this refactor, assuming they are in er_model.er_model and we'll import them.
+from er_model_core.types import Species, ProjectConfig, CarbonConfig
+
+def load_model_config(
+    config_path: Optional[Path] = None, 
+    config_dict: Optional[Dict] = None
+) -> Tuple[List[Species], ProjectConfig, CarbonConfig, str, bool, Dict]:
+    """
+    Loads model configuration from a YAML file or a dictionary.
+
+    Returns:
+        A tuple containing:
+            - List of Species objects
+            - ProjectConfig object
+            - CarbonConfig object
+            - Growth model name (str)
+            - Continuous growth flag (bool)
+            - The raw config dictionary
+    """
+    if config_dict is not None:
+        cfg = config_dict
+    elif config_path is not None:
+        with open(config_path) as f:
+            cfg = yaml.safe_load(f)
+    else:
+        raise ValueError("Either config_path or config_dict must be provided.")
+
+    species_list = [Species(**s) for s in cfg["species"]]
+    project_cfg = ProjectConfig(**cfg["project"])
+    carbon_cfg = CarbonConfig(**cfg["carbon"])
+    
+    growth_model_name = cfg.get('growth_model', 'chapman_richards') 
+    continuous_growth_flag = cfg.get('continuous_growth', False)
+    
+    return species_list, project_cfg, carbon_cfg, growth_model_name, continuous_growth_flag, cfg 

er_model_core/er_model.py

@@ -0,0 +1,348 @@
+"""
+Core implementation of the Emissions Reduction (ER) model for mangrove projects.
+"""
+from dataclasses import dataclass
+from pathlib import Path
+from typing import Dict, List, Optional, Tuple
+
+import numpy as np
+import pandas as pd
+import yaml
+import warnings
+
+from er_model_core.allometry import calculate_biomass
+from er_model_core.metrics import calculate_carbon
+from er_model_core.config_loader import load_model_config
+from er_model_core.types import Species, ProjectConfig, CarbonConfig
+
+# Growth model imports - ensure all are imported
+from er_model_core.growth_models.chapman_richards import chapman_richards_growth
+from er_model_core.growth_models.linear import linear_growth, linear_plateau_growth
+from er_model_core.growth_models.declining_increment import declining_increment_growth, continuous_declining_increment_growth
+
+
+@dataclass
+class Species:
+    """Species-specific parameters for growth and carbon calculations."""
+    name: str
+    planting_density: float
+    # Old style
+    mortality_rates: Optional[Dict[str, float]] = None
+    # New style
+    m_ref: Optional[float] = None
+    DBH_ref: Optional[float] = None
+    p: Optional[float] = None
+    chapman_richards: Dict[str, Dict[str, float]] = None
+    allometry: Dict[str, float] = None
+    initial_values: Dict[str, float] = None
+    linear: Optional[Dict[str, Dict[str, float]]] = None
+    linear_plateau: Optional[Dict[str, Dict[str, float]]] = None
+    declining_increment: Optional[Dict[str, Dict[str, float]]] = None
+
+
+@dataclass
+class ProjectConfig:
+    """Project configuration parameters."""
+    duration_years: int
+    planting_schedule: Dict[str, float]
+
+
+@dataclass
+class CarbonConfig:
+    """Carbon conversion and adjustment parameters."""
+    biomass_to_carbon: float
+    carbon_to_co2: float
+    buffer_percentage: float
+    leakage_percentage: float
+    baseline_emissions: float
+    soil_carbon_per_ha_per_year: float = 0.0
+
+
+class ERModel:
+    """
+    Emissions Reduction Model for mangrove projects.
+    
+    Calculates carbon sequestration over time based on tree growth,
+    mortality, and carbon conversion factors.
+    """
+    
+    def __init__(self, config_path: Optional[Path] = None, config: Optional[dict] = None):
+        """
+        Initialize the model from a YAML configuration file or a config dict.
+        Args:
+            config_path: Path to the YAML configuration file
+            config: Config dict (optional)
+        """
+        (
+            self.species, 
+            self.project, 
+            self.carbon, 
+            self.growth_model, 
+            self.continuous_growth,
+            self.raw_config  # Storing the raw config might be useful for other parts
+        ) = load_model_config(config_path=config_path, config_dict=config)
+        
+        self.results: Optional[pd.DataFrame] = None
+        self.species_results: Optional[pd.DataFrame] = None
+        self.scenario_results: Optional[pd.DataFrame] = None
+        self.scenarios = self.raw_config.get("scenarios", {
+            "area": 1000.0,
+            "dbh_range": [1.0, 20.0],
+            "height_range": [0.5, 12.0],
+            "growth_rate_factor": 1.0
+        })
+        
+    def calculate_cohort_surviving_trees(self, planting_year: int, current_year: int, initial_trees: float, species: Species, plateau_density: Optional[float] = None, growth_model: str = None) -> float:
+        """
+        Calculate surviving trees for a cohort planted in planting_year, in current_year.
+        Uses either per-year mortality (from config) or DBH-dependent mortality.
+        growth_model: which growth model to use for DBH (e.g., 'linear', 'linear_plateau', etc.)
+        """
+        if growth_model is None:
+            growth_model = getattr(self, 'growth_model', 'chapman_richards')
+        age = current_year - planting_year + 1
+        if age < 1:
+            return 0
+        if initial_trees == 0:
+            # Suppress debug output for zero-initial-trees cohorts
+            return 0
+        N_live = initial_trees
+        for year in range(1, age + 1):
+            debug_info = {
+                'planting_year': planting_year,
+                'current_year': current_year,
+                'cohort_age': year,
+                'initial_trees': initial_trees,
+                'N_live_before': N_live
+            }
+            if species.mortality_rates is not None:
+                year_key = f"year_{year}"
+                if year_key in species.mortality_rates:
+                    mort_rate = species.mortality_rates[year_key]
+                else:
+                    mort_rate = species.mortality_rates.get("subsequent", 0)
+                    print(f"[DEBUG] Year key '{year_key}' not found in mortality_rates for {species.name}. Using 'subsequent' or 0. Available keys: {list(species.mortality_rates.keys())}")
+                m = mort_rate / 100.0
+                debug_info['mortality_logic'] = 'annual'
+                debug_info['mortality_rate_percent'] = mort_rate
+            else:
+                growth_func, dbh_params = self.get_growth_function_and_params(species, growth_model, 'dbh')
+                dbh = growth_func(year, dbh_params, species.initial_values["dbh"])
+                dbh = max(dbh, 1.0)
+                m_ref = species.m_ref if species.m_ref is not None else 0.16
+                DBH_ref = species.DBH_ref if species.DBH_ref is not None else 9.0
+                p = species.p if species.p is not None else 1.493
+                m = m_ref * (DBH_ref / dbh) ** p
+                m = min(max(m, 0), 0.99)
+                debug_info['mortality_logic'] = 'dbh-dependent'
+                debug_info['mortality_rate_percent'] = m * 100
+                debug_info['dbh'] = dbh
+            N_live = N_live * (1 - m)
+            debug_info['N_live_after'] = N_live
+            print(f"[DEBUG] {species.name} | PlantingYear: {planting_year} | Year: {current_year} | CohortAge: {year} | MortalityLogic: {debug_info['mortality_logic']} | MortalityRate(%): {debug_info['mortality_rate_percent']} | N_live_before: {debug_info['N_live_before']} | N_live_after: {debug_info['N_live_after']}")
+        return N_live
+
+    def calculate_total_surviving_trees(self, year: int) -> Dict[str, float]:
+        """
+        Calculate total surviving trees for each species in a given year, summing across all cohorts.
+        Returns a dict: {species_name: total_surviving_trees}
+        """
+        growth_model = getattr(self, 'growth_model', 'chapman_richards')
+        totals = {}
+        for species in self.species:
+            total = 0
+            for planting_year, area in self.project.planting_schedule.items():
+                py = int(planting_year.split("_")[1])
+                initial_trees = species.planting_density * area
+                # Use plateau_density as the year-5 value for this cohort
+                plateau_density = species.planting_density * area if 5 <= (year - py + 1) else None
+                total += self.calculate_cohort_surviving_trees(py, year, initial_trees, species, plateau_density, growth_model)
+            totals[species.name] = total
+        return totals
+
+    def calculate_cumulative_area(self, year: int) -> float:
+        """
+        Calculate cumulative area planted up to and including the given year.
+        """
+        total = 0
+        for planting_year, area in self.project.planting_schedule.items():
+            py = int(planting_year.split("_")[1])
+            if py <= year:
+                total += area
+        return total
+
+    def get_growth_function_and_params(self, species: Species, growth_model_name: str, dim: str) -> Tuple[callable, Dict[str, float]]:
+        """
+        Get the growth function and its parameters for a given species, model, and dimension (dbh or height).
+        Uses self.continuous_growth to select continuous version if applicable.
+        """
+        params = species.__dict__.get(growth_model_name, {}).get(dim, {})
+        initial_value_key = f"initial_{dim}"
+        # initial_value = species.initial_values.get(dim, 0) # This was the old way, let's ensure params are complete
+
+        growth_functions = {
+            "chapman_richards": chapman_richards_growth,
+            "linear": linear_growth,
+            "linear_plateau": linear_plateau_growth,
+            "declining_increment": declining_increment_growth,
+            # Potentially add more models here
+        }
+
+        continuous_growth_functions = {
+            "declining_increment": continuous_declining_increment_growth,
+            # Add other continuous versions if they exist
+        }
+
+        if self.continuous_growth and growth_model_name in continuous_growth_functions:
+            selected_func = continuous_growth_functions[growth_model_name]
+        elif growth_model_name in growth_functions:
+            selected_func = growth_functions[growth_model_name]
+        else:
+            # Fallback or error if model name is not found
+            # For now, let's use chapman_richards as a default if self.growth_model itself isn't specific enough
+            # However, growth_model_name parameter should be the specific one.
+            warnings.warn(f"Growth model '{growth_model_name}' not found or not supported. Defaulting to chapman_richards for {species.name} {dim}. Check config.")
+            selected_func = chapman_richards_growth # Default fallback
+            params = species.chapman_richards.get(dim, {}) if species.chapman_richards else {}
+
+        if not params:
+            warnings.warn(f"No parameters found for {growth_model_name} {dim} for species {species.name}. Growth may be incorrect.")
+            # Provide some default empty params if really needed, or let it fail if func expects them.
+
+        return selected_func, params
+
+    def calculate_carbon_for_species(self, species: Species, age: int, area: float, cohort_age: int) -> float:
+        """
+        Calculate carbon sequestration for a single species, cohort, and cohort age.
+        Args:
+            species: Species parameters
+            age: Project year (not used for growth)
+            area: Planted area in hectares
+            cohort_age: Age of this cohort (years since planting)
+        Returns:
+            Carbon sequestration in tCO2
+        """
+        if cohort_age < 1:
+            return 0
+        initial_trees = species.planting_density * area
+        plateau_density = species.planting_density * area if cohort_age >= 5 else None
+        surviving = self.calculate_cohort_surviving_trees(1, cohort_age, initial_trees, species, plateau_density, self.growth_model)
+        dbh_func, dbh_params = self.get_growth_function_and_params(species, self.growth_model, 'dbh')
+        height_func, height_params = self.get_growth_function_and_params(species, self.growth_model, 'height')
+        dbh = dbh_func(cohort_age, dbh_params, species.initial_values["dbh"])
+        height = height_func(cohort_age, height_params, species.initial_values["height"])
+        biomass = calculate_biomass(dbh, height, species.name, species.allometry)
+        carbon = calculate_carbon(
+            biomass * surviving,
+            self.carbon.biomass_to_carbon,
+            self.carbon.carbon_to_co2
+        )
+        return carbon
+
+    def run(self) -> Tuple[pd.DataFrame, pd.DataFrame]:
+        """
+        Execute the full ER calculation pipeline.
+        Returns:
+            Tuple of (yearly results DataFrame, species results DataFrame)
+        """
+        years = range(1, self.project.duration_years + 1)
+        results = []
+        species_results = []
+        species_metrics_rows = []  # For new per-year, per-species metrics
+        cumulative_soil_carbon = 0
+        for year in years:
+            year_results = {"year": year}
+            species_year_results = {"Year": year}
+            total_carbon = 0
+            cumulative_area = self.calculate_cumulative_area(year)
+            for species in self.species:
+                species_carbon = 0
+                # --- New metrics ---
+                total_surviving = 0
+                total_dbh = 0
+                total_height = 0
+                total_biomass_per_tree = 0
+                total_biomass = 0
+                n_cohorts = 0
+                for planting_year, area in self.project.planting_schedule.items():
+                    py = int(planting_year.split("_")[1])
+                    cohort_age = year - py + 1
+                    if cohort_age < 1:
+                        continue
+                    initial_trees = species.planting_density * area
+                    plateau_density = species.planting_density * area if cohort_age >= 5 else None
+                    surviving = self.calculate_cohort_surviving_trees(1, cohort_age, initial_trees, species, plateau_density, self.growth_model)
+                    dbh_func, dbh_params = self.get_growth_function_and_params(species, self.growth_model, 'dbh')
+                    height_func, height_params = self.get_growth_function_and_params(species, self.growth_model, 'height')
+                    dbh = dbh_func(cohort_age, dbh_params, species.initial_values["dbh"])
+                    height = height_func(cohort_age, height_params, species.initial_values["height"])
+                    biomass_per_tree = calculate_biomass(dbh, height, species.name, species.allometry)
+                    total_surviving += surviving
+                    total_dbh += dbh * surviving
+                    total_height += height * surviving
+                    total_biomass_per_tree += biomass_per_tree * surviving
+                    total_biomass += biomass_per_tree * surviving
+                    n_cohorts += surviving
+                    # --- End new metrics ---
+                    # Existing carbon calculation
+                    carbon = self.calculate_carbon_for_species(species, year, area, cohort_age)
+                    species_carbon += carbon
+                total_carbon += species_carbon
+                species_key = f"{species.name} tCO2"
+                species_year_results[species_key] = species_carbon
+                # Store per-year, per-species metrics
+                if total_surviving > 0:
+                    avg_dbh = total_dbh / total_surviving
+                    avg_height = total_height / total_surviving
+                    avg_biomass_per_tree = total_biomass_per_tree / total_surviving
+                else:
+                    avg_dbh = 0
+                    avg_height = 0
+                    avg_biomass_per_tree = 0
+                species_metrics_rows.append({
+                    "Year": year,
+                    "Species": species.name,
+                    "Surviving Trees": total_surviving,
+                    "DBH (cm)": avg_dbh,
+                    "Height (m)": avg_height,
+                    "Biomass per Tree (kg)": avg_biomass_per_tree,
+                    "Total Biomass (kg)": total_biomass
+                })
+            species_year_results["Total tCO2"] = total_carbon
+            species_year_results["Cumulative ha"] = cumulative_area
+            species_year_results["tCO2/ha"] = total_carbon / cumulative_area if cumulative_area > 0 else 0
+            gross_carbon = total_carbon
+            buffer_carbon = gross_carbon * (1 - self.carbon.buffer_percentage / 100)
+            buffer_carbon -= self.carbon.leakage_percentage / 100 * gross_carbon
+            buffer_carbon -= self.carbon.baseline_emissions * cumulative_area
+            # Cumulative soil carbon: add 1 t/ha for every hectare ever planted, each year
+            soil_carbon = 0
+            if hasattr(self.carbon, 'soil_carbon_per_ha_per_year'):
+                soil_carbon = self.carbon.soil_carbon_per_ha_per_year * cumulative_area
+            cumulative_soil_carbon += soil_carbon
+            gross_carbon_with_soil = gross_carbon + cumulative_soil_carbon
+            buffer_carbon_with_soil = buffer_carbon + (cumulative_soil_carbon * (1 - self.carbon.buffer_percentage / 100))
+            year_results.update({
+                "gross_carbon": gross_carbon,
+                "buffer_carbon": buffer_carbon,
+                "cumulative_area": cumulative_area,
+                "gross_carbon_with_soil": gross_carbon_with_soil,
+                "buffer_carbon_with_soil": buffer_carbon_with_soil,
+                "soil_carbon": soil_carbon,
+                "cumulative_soil_carbon": cumulative_soil_carbon
+            })
+            results.append(year_results)
+            species_results.append(species_year_results)
+        self.results = pd.DataFrame(results)
+        self.species_results = pd.DataFrame(species_results)
+        self.species_metrics = pd.DataFrame(species_metrics_rows)
+        return self.results, self.species_results
+
+    def save_results(self, output_path: Path) -> None:
+        """
+        Save the main results DataFrame to a CSV file.
+        """
+        if self.results is None:
+            raise ValueError("Results have not been calculated yet. Run model.run() first.")
+        output_path.parent.mkdir(parents=True, exist_ok=True)
+        self.results.to_csv(output_path, index=False) 

er_model_core/growth_models/chapman_richards.py

@@ -0,0 +1,11 @@
+"""
+Chapman-Richards Growth Model (ARCHIVED/NOT IN USE)
+"""
+import numpy as np
+
+def chapman_richards_growth(age: float, params: dict, initial_value: float) -> float:
+    """
+    Chapman-Richards growth function.
+    """
+    a, b, c = params["a"], params["b"], params["c"]
+    return initial_value + (a - initial_value) * (1 - np.exp(-b * age)) ** c 

er_model_core/growth_models/declining_increment.py

@@ -0,0 +1,31 @@
+"""
+Declining Increment Growth Model (ACTIVE)
+
+Contains both discrete and continuous forms.
+"""
+import numpy as np
+
+def declining_increment_growth(age: float, params: dict, initial_value: float) -> float:
+    """
+    Discrete declining increment growth: sum annual increments, never negative.
+    """
+    r0 = params["r0"]
+    T_m = params["T_m"]
+    total = initial_value
+    for i in range(1, int(np.floor(age)) + 1):
+        increment = r0 * max(0, 1 - (i - 1) / T_m)
+        total += increment
+    frac = age - int(np.floor(age))
+    if frac > 0:
+        i = int(np.floor(age)) + 1
+        increment = r0 * max(0, 1 - (i - 1) / T_m)
+        total += frac * increment
+    return total
+
+def continuous_declining_increment_growth(age: float, params: dict, initial_value: float) -> float:
+    """
+    Continuous declining increment growth: analytical formula.
+    """
+    r0 = params["r0"]
+    T_m = params["T_m"]
+    return initial_value + r0 * (age - age**2 / (2 * T_m)) 

er_model_core/growth_models/linear.py

@@ -0,0 +1,22 @@
+"""
+Linear and Linear Plateau Growth Models (ARCHIVED/NOT IN USE)
+"""
+
+def linear_growth(age: float, params: dict, initial_value: float) -> float:
+    """
+    Simple linear growth function.
+    """
+    r = params["r"]
+    return initial_value + r * age
+
+def linear_plateau_growth(age: float, params: dict, initial_value: float) -> float:
+    """
+    Linear growth with plateau.
+    """
+    r = params["r"]
+    T_p = params["T_p"]
+    a = params["a"]
+    if age <= T_p:
+        return initial_value + r * age
+    else:
+        return initial_value + a 

er_model_core/types.py

@@ -0,0 +1,36 @@
+from dataclasses import dataclass
+from typing import Dict, Optional
+
+@dataclass
+class Species:
+    """Species-specific parameters for growth and carbon calculations."""
+    name: str
+    planting_density: float
+    # Old style
+    mortality_rates: Optional[Dict[str, float]] = None
+    # New style
+    m_ref: Optional[float] = None
+    DBH_ref: Optional[float] = None
+    p: Optional[float] = None
+    chapman_richards: Dict[str, Dict[str, float]] = None
+    allometry: Dict[str, float] = None
+    initial_values: Dict[str, float] = None
+    linear: Optional[Dict[str, Dict[str, float]]] = None
+    linear_plateau: Optional[Dict[str, Dict[str, float]]] = None
+    declining_increment: Optional[Dict[str, Dict[str, float]]] = None
+
+@dataclass
+class ProjectConfig:
+    """Project configuration parameters."""
+    duration_years: int
+    planting_schedule: Dict[str, float]
+
+@dataclass
+class CarbonConfig:
+    """Carbon conversion and adjustment parameters."""
+    biomass_to_carbon: float
+    carbon_to_co2: float
+    buffer_percentage: float
+    leakage_percentage: float
+    baseline_emissions: float
+    soil_carbon_per_ha_per_year: float = 0.0 

models/declining_increment/README.md

@@ -0,0 +1,28 @@
+# Declining Increment Model
+
+This folder documents the **Declining Increment** growth model used in the Mangrove ER Dashboard.
+
+- **Growth Model:** Declining Increment (annual increment declines linearly to zero; total size is the sum of non-negative increments)
+- **Parameters:** See `configs/declining_increment.yaml`
+- **Purpose:** For comparison with other models; may be appropriate for species or sites where growth slows with age.
+- **Usage:** Available as a tab in the dashboard for side-by-side comparison, or as a standalone dashboard.
+
+## Model Details
+
+### Discrete Version (Default)
+- **Logic:** Growth is calculated by summing annual increments for each year of age. The increment for year \(i\) is \( r_0 \cdot \max(0, 1 - \frac{i-1}{T_m}) \).
+- **No Negative Growth:** Increments decline to zero but never go negative, so DBH and height never decrease.
+- **Fractional Years:** If the age is not an integer, a proportional part of the next increment is added.
+- **Mortality:** Annual mortality is compounded each year according to the schedule in the config.
+
+### Continuous Version (Optional)
+- **Logic:** Uses the analytical formula:  
+  \( \text{size}(t) = \text{initial} + r_0 \cdot (t - \frac{t^2}{2T_m}) \)
+- **How to Enable:** Set `continuous_growth: true` in `configs/declining_increment.yaml`.
+- **Implications:** This version is smoother and may yield slightly different results than the discrete sum, especially for short time periods or small Tm.
+
+### Display Names
+- Species are now displayed as "Rhizophora spp." and "Avicennia germinans" throughout the dashboard for clarity.
+
+**Note:**
+This model is flexible and can be toggled between discrete and continuous modes via the config. Small differences in ERs compared to other implementations (e.g., Jupyter notebooks) may arise from differences in how increments, mortality, or fractional years are handled. 

requirements.txt

@@ -0,0 +1,8 @@
+click>=8.1.7
+gradio>=4.0.0
+matplotlib>=3.8.2
+numpy>=1.26.3
+pandas>=2.2.0
+plotly>=5.0.0
+pyyaml>=6.0.1
+scipy>=1.12.0

scripts/run_pipeline.py

@@ -4,14 +4,13 @@ Main pipeline script for running the ER model calculations.
 """
 import sys
 from pathlib import Path
+import yaml
+import click
+from er_model_core.er_model import ERModel
 
 # Add the project root to the Python path
 sys.path.append(str(Path(__file__).parent.parent))
 
-import click
-
-from src.er_model import ERModel
-
 
 @click.command()
 @click.option(

setup.py

@@ -3,7 +3,8 @@ from setuptools import find_packages, setup
 setup(
     name="er_model",
     version="0.1.0",
-    packages=find_packages(),
+    package_dir={'': 'er_model_core'},
+    packages=find_packages(where='er_model_core'),
     install_requires=[
         "pandas>=2.2.0",
         "numpy>=1.26.3",

src/er_model.py

@@ -8,20 +8,32 @@ from typing import Dict, List, Optional, Tuple
 import numpy as np
 import pandas as pd
 import yaml
+import warnings
 
 from .allometry import calculate_biomass
 from .metrics import calculate_carbon
 
+# Growth model imports
+from .growth_models.declining_increment import declining_increment_growth, continuous_declining_increment_growth
+
 
 @dataclass
 class Species:
     """Species-specific parameters for growth and carbon calculations."""
     name: str
     planting_density: float
-    mortality_rates: Dict[str, float]
-    chapman_richards: Dict[str, Dict[str, float]]  # Nested dict for dbh and height parameters
-    allometry: Dict[str, float]
-    initial_values: Dict[str, float]
+    # Old style
+    mortality_rates: Optional[Dict[str, float]] = None
+    # New style
+    m_ref: Optional[float] = None
+    DBH_ref: Optional[float] = None
+    p: Optional[float] = None
+    chapman_richards: Dict[str, Dict[str, float]] = None
+    allometry: Dict[str, float] = None
+    initial_values: Dict[str, float] = None
+    linear: Optional[Dict[str, Dict[str, float]]] = None
+    linear_plateau: Optional[Dict[str, Dict[str, float]]] = None
+    declining_increment: Optional[Dict[str, Dict[str, float]]] = None
 
 
 @dataclass
@@ -39,6 +51,7 @@ class CarbonConfig:
     buffer_percentage: float
     leakage_percentage: float
     baseline_emissions: float
+    soil_carbon_per_ha_per_year: float = 0.0
 
 
 class ERModel:
@@ -49,275 +62,273 @@ class ERModel:
     mortality, and carbon conversion factors.
     """
     
-    def __init__(self, config_path: Path):
+    def __init__(self, config_path: Path = None, config: dict = None):
         """
-        Initialize the model from a YAML configuration file.
-        
+        Initialize the model from a YAML configuration file or a config dict.
         Args:
             config_path: Path to the YAML configuration file
+            config: Config dict (optional)
         """
-        with open(config_path) as f:
-            config = yaml.safe_load(f)
-            
-        self.species = [Species(**s) for s in config["species"]]
-        self.project = ProjectConfig(**config["project"])
-        self.carbon = CarbonConfig(**config["carbon"])
-        
-        # Store scenario parameters if provided
-        self.scenarios = config.get("scenarios", {
-            "area": 1000.0,
-            "dbh_range": [1.0, 20.0],
-            "height_range": [0.5, 12.0],
-            "growth_rate_factor": 1.0
-        })
-        
-        # Initialize results storage
+        if config is not None:
+            cfg = config
+        else:
+            with open(config_path) as f:
+                cfg = yaml.safe_load(f)
+        self.species = [Species(**s) for s in cfg["species"]]
+        self.project = ProjectConfig(**cfg["project"])
+        self.carbon = CarbonConfig(**cfg["carbon"])
         self.results: Optional[pd.DataFrame] = None
         self.species_results: Optional[pd.DataFrame] = None
         self.scenario_results: Optional[pd.DataFrame] = None
+        self.growth_model = cfg.get('growth_model', 'chapman_richards')
+        self.continuous_growth = cfg.get('continuous_growth', False)
         
-    def calculate_surviving_trees(self, initial_trees: float, year: int, 
-                                mortality_rates: Dict[str, float]) -> float:
+    def calculate_cohort_surviving_trees(self, planting_year: int, current_year: int, initial_trees: float, species: Species, plateau_density: Optional[float] = None, growth_model: str = None) -> float:
         """
-        Calculate surviving trees for a given year based on mortality rates.
-        
-        Args:
-            initial_trees: Initial number of trees planted
-            year: Years since planting (1-based)
-            mortality_rates: Dictionary of mortality rates by year
-            
-        Returns:
-            Number of surviving trees
+        Calculate surviving trees for a cohort planted in planting_year, in current_year.
+        Uses either per-year mortality (from config) or DBH-dependent mortality.
+        growth_model: which growth model to use for DBH (e.g., 'linear', 'linear_plateau', etc.)
         """
-        surviving = initial_trees
-        
-        for y in range(1, year + 1):
-            rate = mortality_rates.get(f"year_{y}", 
-                                    mortality_rates["subsequent"]) / 100.0
-            surviving *= (1 - rate)
-            
-        return surviving
-    
-    def chapman_richards_growth(self, age: float, params: Dict[str, float], initial_value: float) -> float:
+        if growth_model is None:
+            growth_model = getattr(self, 'growth_model', 'chapman_richards')
+        age = current_year - planting_year + 1
+        if age < 1:
+            return 0
+        if initial_trees == 0:
+            # Suppress debug output for zero-initial-trees cohorts
+            return 0
+        N_live = initial_trees
+        for year in range(1, age + 1):
+            debug_info = {
+                'planting_year': planting_year,
+                'current_year': current_year,
+                'cohort_age': year,
+                'initial_trees': initial_trees,
+                'N_live_before': N_live
+            }
+            if species.mortality_rates is not None:
+                year_key = f"year_{year}"
+                if year_key in species.mortality_rates:
+                    mort_rate = species.mortality_rates[year_key]
+                else:
+                    mort_rate = species.mortality_rates.get("subsequent", 0)
+                    print(f"[DEBUG] Year key '{year_key}' not found in mortality_rates for {species.name}. Using 'subsequent' or 0. Available keys: {list(species.mortality_rates.keys())}")
+                m = mort_rate / 100.0
+                debug_info['mortality_logic'] = 'annual'
+                debug_info['mortality_rate_percent'] = mort_rate
+            else:
+                growth_func, dbh_params = self.get_growth_function_and_params(species, growth_model, 'dbh')
+                dbh = growth_func(year, dbh_params, species.initial_values["dbh"])
+                dbh = max(dbh, 1.0)
+                m_ref = species.m_ref if species.m_ref is not None else 0.16
+                DBH_ref = species.DBH_ref if species.DBH_ref is not None else 9.0
+                p = species.p if species.p is not None else 1.493
+                m = m_ref * (DBH_ref / dbh) ** p
+                m = min(max(m, 0), 0.99)
+                debug_info['mortality_logic'] = 'dbh-dependent'
+                debug_info['mortality_rate_percent'] = m * 100
+                debug_info['dbh'] = dbh
+            N_live = N_live * (1 - m)
+            debug_info['N_live_after'] = N_live
+            print(f"[DEBUG] {species.name} | PlantingYear: {planting_year} | Year: {current_year} | CohortAge: {year} | MortalityLogic: {debug_info['mortality_logic']} | MortalityRate(%): {debug_info['mortality_rate_percent']} | N_live_before: {debug_info['N_live_before']} | N_live_after: {debug_info['N_live_after']}")
+        return N_live
+
+    def calculate_total_surviving_trees(self, year: int) -> Dict[str, float]:
         """
-        Calculate growth using Chapman-Richards growth equation.
-        
-        Args:
-            age: Tree age in years
-            params: Dictionary with a, b, c parameters
-            initial_value: Initial value (DBH or height)
-            
-        Returns:
-            Current size (DBH in cm or height in m)
+        Calculate total surviving trees for each species in a given year, summing across all cohorts.
+        Returns a dict: {species_name: total_surviving_trees}
         """
+        growth_model = getattr(self, 'growth_model', 'chapman_richards')
+        totals = {}
+        for species in self.species:
+            total = 0
+            for planting_year, area in self.project.planting_schedule.items():
+                py = int(planting_year.split("_")[1])
+                initial_trees = species.planting_density * area
+                # Use plateau_density as the year-5 value for this cohort
+                plateau_density = species.planting_density * area if 5 <= (year - py + 1) else None
+                total += self.calculate_cohort_surviving_trees(py, year, initial_trees, species, plateau_density, growth_model)
+            totals[species.name] = total
+        return totals
+
+    def calculate_cumulative_area(self, year: int) -> float:
+        """
+        Calculate cumulative area planted up to and including the given year.
+        """
+        total = 0
+        for planting_year, area in self.project.planting_schedule.items():
+            py = int(planting_year.split("_")[1])
+            if py <= year:
+                total += area
+        return total
+
+    @staticmethod
+    def chapman_richards_growth(age: float, params: Dict[str, float], initial_value: float) -> float:
         a, b, c = params["a"], params["b"], params["c"]
         return initial_value + (a - initial_value) * (1 - np.exp(-b * age)) ** c
 
-    def calculate_carbon_for_species(self, species: Species, age: int, area: float) -> float:
+    @staticmethod
+    def linear_growth(age: float, params: Dict[str, float], initial_value: float) -> float:
+        r = params["r"]
+        return initial_value + r * age
+
+    @staticmethod
+    def linear_plateau_growth(age: float, params: Dict[str, float], initial_value: float) -> float:
+        r = params["r"]
+        T_p = params["T_p"]
+        a = params["a"]
+        if age <= T_p:
+            return initial_value + r * age
+        else:
+            return initial_value + a
+
+    @staticmethod
+    def declining_increment_growth(age: float, params: Dict[str, float], initial_value: float) -> float:
+        r0 = params["r0"]
+        T_m = params["T_m"]
+        # Accumulate annual increments, never negative
+        total = initial_value
+        for i in range(1, int(np.floor(age)) + 1):
+            increment = r0 * max(0, 1 - (i - 1) / T_m)
+            total += increment
+        # If age is fractional, add partial increment for the last year
+        frac = age - int(np.floor(age))
+        if frac > 0:
+            i = int(np.floor(age)) + 1
+            increment = r0 * max(0, 1 - (i - 1) / T_m)
+            total += frac * increment
+        return total
+
+    @staticmethod
+    def continuous_declining_increment_growth(age: float, params: Dict[str, float], initial_value: float) -> float:
+        r0 = params["r0"]
+        T_m = params["T_m"]
+        # Continuous formula: initial + r0 * (age - age^2/(2*Tm))
+        return initial_value + r0 * (age - age**2 / (2 * T_m))
+
+    def calculate_carbon_for_species(self, species: Species, age: int, area: float, cohort_age: int) -> float:
         """
-        Calculate carbon sequestration for a single species and age.
-        
+        Calculate carbon sequestration for a single species, cohort, and cohort age.
         Args:
             species: Species parameters
-            age: Age of trees in years
+            age: Project year (not used for growth)
             area: Planted area in hectares
-            
+            cohort_age: Age of this cohort (years since planting)
         Returns:
             Carbon sequestration in tCO2
         """
-        # Calculate surviving trees
+        if cohort_age < 1:
+            return 0
         initial_trees = species.planting_density * area
-        surviving = self.calculate_surviving_trees(
-            initial_trees, age, species.mortality_rates
-        )
-        
-        # Calculate DBH and height using Chapman-Richards
-        dbh = self.chapman_richards_growth(
-            age, 
-            species.chapman_richards["dbh"],
-            species.initial_values["dbh"]
-        )
-        height = self.chapman_richards_growth(
-            age,
-            species.chapman_richards["height"],
-            species.initial_values["height"]
-        )
-        
-        # Calculate biomass using both DBH and height
+        plateau_density = species.planting_density * area if cohort_age >= 5 else None
+        surviving = self.calculate_cohort_surviving_trees(1, cohort_age, initial_trees, species, plateau_density, self.growth_model)
+        dbh_func, dbh_params = self.get_growth_function_and_params(species, self.growth_model, 'dbh')
+        height_func, height_params = self.get_growth_function_and_params(species, self.growth_model, 'height')
+        dbh = dbh_func(cohort_age, dbh_params, species.initial_values["dbh"])
+        height = height_func(cohort_age, height_params, species.initial_values["height"])
         biomass = calculate_biomass(dbh, height, species.name, species.allometry)
-        
-        # Convert to carbon
         carbon = calculate_carbon(
             biomass * surviving,
             self.carbon.biomass_to_carbon,
             self.carbon.carbon_to_co2
         )
-        
         return carbon
-    
-    def run_scenario(self, area: float, dbh_range: List[float], height_range: List[float], 
-                    growth_rate_factor: float) -> pd.Series:
-        """
-        Run a scenario with modified parameters.
-        
-        Args:
-            area: Area to plant in hectares
-            dbh_range: [min_dbh, max_dbh] in cm
-            height_range: [min_height, max_height] in m
-            growth_rate_factor: Factor to multiply growth rate by
-            
-        Returns:
-            Series with carbon sequestration at milestone years
-        """
-        # Create modified species for scenario
-        scenario_species = []
-        for sp in self.species:
-            mod_sp = Species(
-                name=sp.name,
-                planting_density=sp.planting_density,
-                mortality_rates=sp.mortality_rates.copy(),
-                chapman_richards={
-                    "dbh": {
-                        "a": np.mean(dbh_range),
-                        "b": sp.chapman_richards["dbh"]["b"] * growth_rate_factor,
-                        "c": sp.chapman_richards["dbh"]["c"]
-                    },
-                    "height": {
-                        "a": np.mean(height_range),
-                        "b": sp.chapman_richards["height"]["b"] * growth_rate_factor,
-                        "c": sp.chapman_richards["height"]["c"]
-                    }
-                },
-                allometry=sp.allometry.copy(),
-                initial_values=sp.initial_values.copy()
-            )
-            scenario_species.append(mod_sp)
-        
-        # Calculate carbon at milestone years
-        milestones = [5, 10, 15, 30]
-        carbon_values = []
-        
-        for year in milestones:
-            total_carbon = sum(
-                self.calculate_carbon_for_species(sp, year, area)
-                for sp in scenario_species
-            )
-            carbon_values.append(total_carbon)
-        
-        return pd.Series(
-            carbon_values,
-            index=[f"Year {y} tCO2" for y in milestones]
-        )
-    
-    def run(self) -> Tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame]:
+
+    def run(self) -> Tuple[pd.DataFrame, pd.DataFrame]:
         """
         Execute the full ER calculation pipeline.
-        
         Returns:
-            Tuple of (yearly results DataFrame, species results DataFrame, scenario results DataFrame)
+            Tuple of (yearly results DataFrame, species results DataFrame)
         """
         years = range(1, self.project.duration_years + 1)
         results = []
         species_results = []
-        
+        species_metrics_rows = []  # For new per-year, per-species metrics
         for year in years:
             year_results = {"year": year}
             species_year_results = {"Year": year}
             total_carbon = 0
-            cumulative_area = 0
-            
-            # Calculate for each planting cohort and species
-            for planting_year, area in self.project.planting_schedule.items():
-                py = int(planting_year.split("_")[1])
-                if py > year:
-                    continue
-                
-                cumulative_area += area
-                age = year - py + 1
-                
-                for species in self.species:
-                    carbon = self.calculate_carbon_for_species(species, age, area)
-                    total_carbon += carbon
-                    
-                    # Track per-species carbon
-                    species_key = f"{species.name} tCO2"
-                    species_year_results[species_key] = species_year_results.get(species_key, 0) + carbon
-            
-            # Add total carbon and area metrics
+            cumulative_area = self.calculate_cumulative_area(year)
+            for species in self.species:
+                species_carbon = 0
+                # --- New metrics ---
+                total_surviving = 0
+                total_dbh = 0
+                total_height = 0
+                total_biomass_per_tree = 0
+                total_biomass = 0
+                n_cohorts = 0
+                for planting_year, area in self.project.planting_schedule.items():
+                    py = int(planting_year.split("_")[1])
+                    cohort_age = year - py + 1
+                    if cohort_age < 1:
+                        continue
+                    initial_trees = species.planting_density * area
+                    plateau_density = species.planting_density * area if cohort_age >= 5 else None
+                    surviving = self.calculate_cohort_surviving_trees(1, cohort_age, initial_trees, species, plateau_density, self.growth_model)
+                    dbh_func, dbh_params = self.get_growth_function_and_params(species, self.growth_model, 'dbh')
+                    height_func, height_params = self.get_growth_function_and_params(species, self.growth_model, 'height')
+                    dbh = dbh_func(cohort_age, dbh_params, species.initial_values["dbh"])
+                    height = height_func(cohort_age, height_params, species.initial_values["height"])
+                    biomass_per_tree = calculate_biomass(dbh, height, species.name, species.allometry)
+                    total_surviving += surviving
+                    total_dbh += dbh * surviving
+                    total_height += height * surviving
+                    total_biomass_per_tree += biomass_per_tree * surviving
+                    total_biomass += biomass_per_tree * surviving
+                    n_cohorts += surviving
+                    # --- End new metrics ---
+                    # Existing carbon calculation
+                    carbon = self.calculate_carbon_for_species(species, year, area, cohort_age)
+                    species_carbon += carbon
+                total_carbon += species_carbon
+                species_key = f"{species.name} tCO2"
+                species_year_results[species_key] = species_carbon
+                # Store per-year, per-species metrics
+                if total_surviving > 0:
+                    avg_dbh = total_dbh / total_surviving
+                    avg_height = total_height / total_surviving
+                    avg_biomass_per_tree = total_biomass_per_tree / total_surviving
+                else:
+                    avg_dbh = 0
+                    avg_height = 0
+                    avg_biomass_per_tree = 0
+                species_metrics_rows.append({
+                    "Year": year,
+                    "Species": species.name,
+                    "Surviving Trees": total_surviving,
+                    "DBH (cm)": avg_dbh,
+                    "Height (m)": avg_height,
+                    "Biomass per Tree (kg)": avg_biomass_per_tree,
+                    "Total Biomass (kg)": total_biomass
+                })
             species_year_results["Total tCO2"] = total_carbon
             species_year_results["Cumulative ha"] = cumulative_area
             species_year_results["tCO2/ha"] = total_carbon / cumulative_area if cumulative_area > 0 else 0
-            
-            # Apply adjustments for main results
             gross_carbon = total_carbon
-            net_carbon = gross_carbon * (1 - self.carbon.buffer_percentage / 100)
-            net_carbon -= self.carbon.leakage_percentage / 100 * gross_carbon
-            net_carbon -= self.carbon.baseline_emissions * cumulative_area
-            
+            buffer_carbon = gross_carbon * (1 - self.carbon.buffer_percentage / 100)
+            buffer_carbon -= self.carbon.leakage_percentage / 100 * gross_carbon
+            buffer_carbon -= self.carbon.baseline_emissions * cumulative_area
+            # Cumulative soil carbon: add 1 t/ha for every hectare ever planted, each year
+            soil_carbon = 0
+            if hasattr(self.carbon, 'soil_carbon_per_ha_per_year'):
+                soil_carbon = self.carbon.soil_carbon_per_ha_per_year * cumulative_area
+            gross_carbon_with_soil = gross_carbon + soil_carbon
+            buffer_carbon_with_soil = buffer_carbon + soil_carbon
             year_results.update({
                 "gross_carbon": gross_carbon,
-                "net_carbon": net_carbon,
-                f"area_year_{year}": area if f"year_{year}" in self.project.planting_schedule else 0
+                "buffer_carbon": buffer_carbon,
+                "cumulative_area": cumulative_area,
+                "gross_carbon_with_soil": gross_carbon_with_soil,
+                "buffer_carbon_with_soil": buffer_carbon_with_soil,
+                "soil_carbon": soil_carbon
             })
-            
             results.append(year_results)
             species_results.append(species_year_results)
-        
         self.results = pd.DataFrame(results)
         self.species_results = pd.DataFrame(species_results)
-        
-        # Run scenarios
-        scenarios = []
-        dbh_min, dbh_max = self.scenarios["dbh_range"]
-        height_min, height_max = self.scenarios["height_range"]
-        area = self.scenarios["area"]
-        growth_factor = self.scenarios["growth_rate_factor"]
-        
-        # Base scenario
-        base_scenario = self.run_scenario(
-            area, [dbh_min, dbh_max], [height_min, height_max], growth_factor
-        )
-        scenarios.append({
-            "Scenario": "Base",
-            "DBH Range": f"{dbh_min:.1f}-{dbh_max:.1f}",
-            "Height Range": f"{height_min:.1f}-{height_max:.1f}",
-            "Growth Rate": f"{growth_factor:.1f}x",
-            **base_scenario
-        })
-        
-        # Low growth scenario
-        low_scenario = self.run_scenario(
-            area, 
-            [dbh_min * 0.5, dbh_max * 0.5],
-            [height_min * 0.5, height_max * 0.5],
-            growth_factor
-        )
-        scenarios.append({
-            "Scenario": "Low Growth",
-            "DBH Range": f"{dbh_min*0.5:.1f}-{dbh_max*0.5:.1f}",
-            "Height Range": f"{height_min*0.5:.1f}-{height_max*0.5:.1f}",
-            "Growth Rate": f"{growth_factor:.1f}x",
-            **low_scenario
-        })
-        
-        # High growth scenario
-        high_scenario = self.run_scenario(
-            area,
-            [dbh_min * 1.5, dbh_max * 1.5],
-            [height_min * 1.5, height_max * 1.5],
-            growth_factor
-        )
-        scenarios.append({
-            "Scenario": "High Growth",
-            "DBH Range": f"{dbh_min*1.5:.1f}-{dbh_max*1.5:.1f}",
-            "Height Range": f"{height_min*1.5:.1f}-{height_max*1.5:.1f}",
-            "Growth Rate": f"{growth_factor:.1f}x",
-            **high_scenario
-        })
-        
-        self.scenario_results = pd.DataFrame(scenarios)
-        
-        return self.results, self.species_results, self.scenario_results
+        self.species_metrics = pd.DataFrame(species_metrics_rows)
+        return self.results, self.species_results
 
     def save_results(self, output_path: Path) -> None:
         """
@@ -328,4 +339,62 @@ class ERModel:
         """
         if self.results is None:
             raise ValueError("No results available. Run the model first.")
-        self.results.to_csv(output_path, index=False) 
+        self.results.to_csv(output_path, index=False)
+
+    def get_growth_function_and_params(self, species, growth_model, dim):
+        """
+        Returns the correct growth function and parameter dict for the given species and dimension (dbh or height).
+        """
+        if growth_model == "linear":
+            # from growth_models.linear import linear_growth
+            func = None  # ARCHIVED/NOT IN USE
+            params = species.linear[dim]
+        elif growth_model == "linear_plateau":
+            # from growth_models.linear import linear_plateau_growth
+            func = None  # ARCHIVED/NOT IN USE
+            params = species.linear_plateau[dim]
+        elif growth_model == "declining_increment":
+            if getattr(self, 'continuous_growth', False):
+                func = continuous_declining_increment_growth
+            else:
+                func = declining_increment_growth
+            params = species.declining_increment[dim]
+        else:
+            # from growth_models.chapman_richards import chapman_richards_growth
+            func = None  # ARCHIVED/NOT IN USE
+            params = species.chapman_richards[dim]
+        return func, params
+
+# --- Parameter sweep/test for plausible survival curves ---
+def test_dbh_mortality_sweep():
+    import matplotlib.pyplot as plt
+    import numpy as np
+    m_refs = [0.01, 0.05, 0.1, 0.16]
+    ps = [1.0, 1.5, 2.0]
+    DBH_ref = 9.0
+    years = np.arange(1, 31)
+    initial_trees = 1000
+    results = {}
+    for m_ref in m_refs:
+        for p in ps:
+            N_live = initial_trees
+            N_lives = []
+            for year in years:
+                dbh = 1.0 + (year - 1) * 0.5  # simple linear DBH growth for test
+                dbh = max(dbh, 1.0)
+                m = m_ref * (DBH_ref / dbh) ** p
+                m = min(max(m, 0), 0.99)
+                N_live = N_live * (1 - m)
+                N_lives.append(N_live)
+            results[(m_ref, p)] = N_lives
+    plt.figure(figsize=(10,6))
+    for (m_ref, p), N_lives in results.items():
+        plt.plot(years, N_lives, label=f"m_ref={m_ref}, p={p}")
+    plt.xlabel("Year")
+    plt.ylabel("Surviving Trees")
+    plt.title("DBH-dependent Mortality Parameter Sweep")
+    plt.legend()
+    plt.grid(True)
+    plt.show()
+
+# To run the test, call test_dbh_mortality_sweep() from __main__ or a notebook. 

tests/test_er_model.py

@@ -6,10 +6,11 @@ from pathlib import Path
 import numpy as np
 import pytest
 import yaml
+import pandas as pd
 
-from src.allometry import calculate_biomass
-from src.er_model import ERModel, Species
-from src.metrics import calculate_carbon
+from er_model_core.allometry import calculate_biomass
+from er_model_core.er_model import ERModel, Species
+from er_model_core.metrics import calculate_carbon
 
 
 @pytest.fixture
@@ -20,20 +21,18 @@ def sample_config():
             {
                 "name": "Test Species",
                 "planting_density": 1000,
-                "mortality_rates": {
-                    "year_1": 10,
-                    "year_2": 5,
-                    "subsequent": 2
-                },
+                "m_ref": 0.16,
+                "DBH_ref": 9.0,
+                "p": 1.493,
                 "chapman_richards": {
-                    "a": 30,
-                    "b": 0.15,
-                    "c": 1.5
+                    "dbh": {"a": 30, "b": 0.15, "c": 1.5},
+                    "height": {"a": 10, "b": 0.05, "c": 1.2}
                 },
                 "allometry": {
                     "biomass_equation": "0.2 * (dbh ** 2.4)",
                     "root_shoot_ratio": 0.4
-                }
+                },
+                "initial_values": {"dbh": 1.0, "height": 0.5}
             }
         ],
         "project": {
@@ -71,18 +70,16 @@ def test_model_initialization(config_file):
 
 
 def test_surviving_trees_calculation(config_file):
-    """Test tree survival calculations."""
+    """Test tree survival calculations with DBH-dependent mortality."""
     model = ERModel(config_file)
     species = model.species[0]
-    
-    # Test first year mortality
-    surviving = model.calculate_surviving_trees(1000, 1, species.mortality_rates)
-    assert surviving == 900  # 1000 * (1 - 0.10)
-    
-    # Test cumulative mortality
-    surviving = model.calculate_surviving_trees(1000, 2, species.mortality_rates)
-    expected = 1000 * (1 - 0.10) * (1 - 0.05)
-    assert np.isclose(surviving, expected)
+    # Test first year mortality (should use DBH-dependent formula)
+    dbh = 9.0
+    m = species.m_ref * (species.DBH_ref / dbh) ** species.p
+    initial_trees = 1000
+    N_live = initial_trees * (1 - m)
+    expected = initial_trees * (1 - 0.16)
+    assert np.isclose(N_live, expected, atol=1e-3)
 
 
 def test_chapman_richards_growth(config_file):
@@ -132,4 +129,112 @@ def test_full_pipeline(config_file):
     assert "year" in results.columns
     assert "gross_carbon" in results.columns
     assert "net_carbon" in results.columns
-    assert all(results["net_carbon"] <= results["gross_carbon"]) 
+    assert all(results["net_carbon"] <= results["gross_carbon"])
+
+
+def test_dbh_dependent_mortality():
+    """Test that DBH-dependent mortality matches expected at DBH=9 and year-5 plateau."""
+    from er_model_core.er_model import ERModel, Species
+    # Minimal species config
+    species = Species(
+        name="Test Species",
+        planting_density=1000,
+        m_ref=0.16,
+        DBH_ref=9.0,
+        p=1.493,
+        chapman_richards={
+            "dbh": {"a": 9.0, "b": 0.5, "c": 1.0},
+            "height": {"a": 5.0, "b": 0.2, "c": 1.0}
+        },
+        allometry={"equation": "0.2 * (dbh ** 2.4)", "root_shoot_ratio": 0.4},
+        initial_values={"dbh": 9.0, "height": 5.0}
+    )
+    # At DBH = 9, mortality = 0.16
+    dbh = 9.0
+    m = species.m_ref * (species.DBH_ref / dbh) ** species.p
+    assert np.isclose(m, 0.16, atol=1e-3)
+
+    # Simulate 5 years, enforce plateau at year 5
+    initial_trees = 1000
+    plateau_density = 2000
+    class DummyModel:
+        def chapman_richards_growth(self, age, params, initial):
+            return dbh  # Always 9 for this test
+    model = DummyModel()
+    N_live = initial_trees
+    for y in range(1, 6):
+        m = species.m_ref * (species.DBH_ref / dbh) ** species.p
+        N_live = N_live * (1 - m)
+        if y == 5:
+            N_live = plateau_density
+    assert N_live == plateau_density
+
+
+def test_allometric_equation_species_A_B():
+    """Test Zanvo et al. 2023 allometric equations for both species."""
+    from er_model_core.allometry import calculate_biomass
+    # Test values
+    dbh = 10.0  # cm
+    height = 5.0  # m
+    # species_A (Rhizophora spp.)
+    expected_A = 2.0738 * (dbh**2 * height)**0.67628
+    result_A = calculate_biomass(dbh, height, "species_A", {})
+    assert np.isclose(result_A, expected_A, rtol=1e-6), f"species_A: got {result_A}, expected {expected_A}"
+    # species_B (Avicennia germinans)
+    expected_B = 1.5595 * (dbh**2 * height)**0.55864
+    result_B = calculate_biomass(dbh, height, "species_B", {})
+    assert np.isclose(result_B, expected_B, rtol=1e-6), f"species_B: got {result_B}, expected {expected_B}" 
+
+# --- Parameter sweep/test for plausible survival curves (moved from er_model.py) ---
+def test_dbh_mortality_sweep():
+    import matplotlib.pyplot as plt
+    import numpy as np
+    m_refs = [0.01, 0.05, 0.1, 0.16]
+    ps = [1.0, 1.5, 2.0]
+    DBH_ref = 9.0
+    years = np.arange(1, 31)
+    initial_trees = 1000
+    results = {}
+    for m_ref in m_refs:
+        for p in ps:
+            N_live = initial_trees
+            N_lives = []
+            for year in years:
+                # This test requires a simplified DBH growth for its internal logic.
+                # It does not use the main ERModel's growth functions directly.
+                dbh_test_growth = 1.0 + (year - 1) * 0.5  # simple linear DBH growth for test
+                dbh = max(dbh_test_growth, 1.0)
+                
+                # Mortality calculation as per ERModel logic being tested
+                m = m_ref * (DBH_ref / dbh) ** p
+                m = min(max(m, 0), 0.99) 
+                N_live = N_live * (1 - m)
+                N_lives.append(N_live)
+            results[(m_ref, p)] = N_lives
+
+    # This test originally had plotting. For automated tests, assertions are better.
+    # If visual inspection is needed, this part can be run in a notebook.
+    # For now, let's assert that the number of trees is non-increasing.
+    for (m_ref, p), N_lives_curve in results.items():
+        for i in range(len(N_lives_curve) - 1):
+            assert N_lives_curve[i+1] <= N_lives_curve[i], f"Mortality sweep {m_ref, p}: Trees increased from {N_lives_curve[i]} to {N_lives_curve[i+1]}"
+        assert N_lives_curve[-1] < initial_trees, f"Mortality sweep {m_ref, p}: No mortality occurred or trees increased."
+
+    # plt.figure(figsize=(10,6))
+    # for (m_ref, p), N_lives_plot in results.items():
+    #     plt.plot(years, N_lives_plot, label=f"m_ref={m_ref}, p={p}")
+    # plt.xlabel("Year")
+    # plt.ylabel("Surviving Trees")
+    # plt.title("DBH-dependent Mortality Parameter Sweep")
+    # plt.legend()
+    # plt.grid(True)
+    # To save plot instead of showing:
+    # output_dir = Path(__file__).parent.parent / "outputs" / "test_plots"
+    # output_dir.mkdir(parents=True, exist_ok=True)
+    # plt.savefig(output_dir / "dbh_mortality_sweep.png")
+    # plt.close()
+    print("test_dbh_mortality_sweep completed basic assertions.")
+
+# If you want to run this specific test and see the plot (e.g., during development):
+# if __name__ == "__main__":
+#     test_dbh_mortality_sweep() # You would need to handle imports for ERModel parts if called directly here