Gittyup keys
# ============================================================
# MONAD-FIELD SPARC – DWARF GALAXY FIT (DDO154)
# Tests universality of γ in low‑baryon regime.
# ============================================================
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from scipy.interpolate import interp1d
# ------------------------------------------------------------
# Load your uploaded SPARC file (already done in previous code)
# Instead of re-loading, assume df exists. If not, use the code from earlier.
# ------------------------------------------------------------
# ------------------------------------------------------------
# 2. Select galaxy (change GALAXY name as needed)
# ------------------------------------------------------------
GALAXY = "DDO154" # classic dark‑matter‑dominated dwarf
# Alternative: "DDO52", "NGC2915", "IC2574", "UGC4325"
df_gal = df[df['ID'].str.strip() == GALAXY].copy()
if len(df_gal) == 0:
raise ValueError(f"Galaxy {GALAXY} not found. Available: {df['ID'].unique()[:20]}...")
print(f"\nLoaded {len(df_gal)} data points for {GALAXY}")
# ------------------------------------------------------------
# 3. Compute baryonic and residual rotation curves (same as before)
# ------------------------------------------------------------
df_gal['Vbary_sq'] = df_gal['Vgas']**2 + df_gal['Vdisk']**2 + df_gal['Vbul']**2
df_gal['Vbary'] = np.sqrt(df_gal['Vbary_sq'])
df_gal['Vobs_sq'] = df_gal['Vobs']**2
df_gal['Vsub_sq'] = df_gal['Vobs_sq'] - df_gal['Vbary_sq']
df_gal['Vsub_sq'] = df_gal['Vsub_sq'].clip(lower=0)
df_gal['Vsub'] = np.sqrt(df_gal['Vsub_sq'])
df_gal['e_Vsub'] = df_gal['Vobs'] * df_gal['e_Vobs'] / df_gal['Vsub']
df_gal['e_Vsub'] = df_gal['e_Vsub'].replace([np.inf, -np.inf], 0).fillna(0)
# ------------------------------------------------------------
# 4. Fit Monad‑Field static‑limit model
# ------------------------------------------------------------
def sub_powerlaw(R, V0, gamma, R0=1.0):
return V0 * (R / R0)**gamma
mask = (df_gal['R'] > 0) & (df_gal['Vsub'] > 0) & (df_gal['e_Vsub'] > 0)
R_fit = df_gal.loc[mask, 'R'].values
Vsub_fit = df_gal.loc[mask, 'Vsub'].values
e_fit = df_gal.loc[mask, 'e_Vsub'].values
try:
popt, pcov = curve_fit(sub_powerlaw, R_fit, Vsub_fit,
p0=[20.0, 0.5], sigma=e_fit, absolute_sigma=True,
bounds=([0, -0.5], [200, 2.0]))
V0_fit, gamma_fit = popt
perr = np.sqrt(np.diag(pcov))
print(f"\nFitted Monad‑Field substrate parameters for {GALAXY}:")
print(f" V0 = {V0_fit:.2f} ± {perr[0]:.2f} km/s")
print(f" gamma = {gamma_fit:.3f} ± {perr[1]:.3f}")
except Exception as e:
print(f"Fit failed: {e}")
V0_fit, gamma_fit = np.nan, np.nan
# ------------------------------------------------------------
# 5. Plot
# ------------------------------------------------------------
plt.figure(figsize=(10, 7))
plt.errorbar(df_gal['R'], df_gal['Vobs'], yerr=df_gal['e_Vobs'],
fmt='o', color='black', capsize=2, label='Observed')
plt.plot(df_gal['R'], df_gal['Vbary'], 'b-', lw=2, label='Baryonic (Newtonian)')
if not np.isnan(V0_fit):
R_model = np.linspace(df_gal['R'].min(), df_gal['R'].max(), 200)
Vsub_model = sub_powerlaw(R_model, V0_fit, gamma_fit)
f_bary = interp1d(df_gal['R'], df_gal['Vbary'], kind='cubic', fill_value='extrapolate')
Vbary_model = f_bary(R_model)
Vtotal_model = np.sqrt(Vbary_model**2 + Vsub_model**2)
plt.plot(R_model, Vtotal_model, 'r--', lw=2,
label=f'Monad‑Field total (γ = {gamma_fit:.2f})')
plt.xlabel('Radius (kpc)')
plt.ylabel('Circular velocity (km/s)')
plt.title(f'Rotation curve of {GALAXY}')
plt.legend()
plt.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig(f'MonadField_{GALAXY}_fit.png', dpi=150)
plt.show()
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