import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split

# -------------------------------
# Example training data (synthetic)
# Replace with experimental dataset: [(frequency, nodal_count)]
# -------------------------------
frequencies = np.linspace(100, 1000, 100)
nodal_counts = 3*np.sin(0.01*frequencies) + 0.5*np.cos(0.005*frequencies) + np.random.normal(0,0.05,len(frequencies))

X = np.vstack([
    frequencies,
    frequencies**2,
    np.sin(frequencies/100),
    np.cos(frequencies/200),
    np.power(1.618, frequencies/500),  # φ-scaling feature
    np.sqrt(frequencies),
    np.log1p(frequencies)
]).T

Y = np.vstack([
    nodal_counts + 0.1*np.random.randn(len(frequencies)),
    nodal_counts*2 + 0.2*np.random.randn(len(frequencies)),
    nodal_counts*0.7 + 0.1*np.random.randn(len(frequencies)),
    nodal_counts*1.1 + 0.2*np.random.randn(len(frequencies))
]).T

# Train/test split
X_train, X_test, Y_train, Y_test = train_test_split(X, Y, test_size=0.2, random_state=42)

# Linear regression fit
model = LinearRegression().fit(X_train, Y_train)

print("Learned T matrix:\n", model.coef_)
print("Intercept:\n", model.intercept_)
print("R^2 test score:", model.score(X_test, Y_test))

# -------------------------------
# φ-Hash conditional tuning functions
# -------------------------------
def phi_hash(x, phi=1.6180339887):
    """Nonlinear φ-hash of feature vector x."""
    return np.mod(phi ** x, 1.0)

def conditional_phi_tune(y_pred, x, y_true, lam=0.01):
    """Apply φ-hashing only if it improves prediction."""
    y_candidate = y_pred + lam * phi_hash(x)
    error_orig = np.linalg.norm(y_pred - y_true)
    error_candidate = np.linalg.norm(y_candidate - y_true)
    if error_candidate < error_orig:
        return y_candidate
    else:
        return y_pred

# Predict parameters at a given frequency
f_query = 432  # example sacred tone
features = np.array([
    f_query,
    f_query**2,
    np.sin(f_query/100),
    np.cos(f_query/200),
    np.power(1.618, f_query/500),
    np.sqrt(f_query),
    np.log1p(f_query)
]).reshape(1,-1)

# Initial prediction from linear model
params = model.predict(features)[0]

# Optional: tune with φ-hashing using test set as reference
# For demonstration, using Y_test mean as pseudo ground truth
y_true_ref = Y_test.mean(axis=0)
params_tuned = conditional_phi_tune(params, features.flatten(), y_true_ref)

cymatic_params = {
    "alpha": abs(params_tuned[0]),
    "beta": abs(params_tuned[1]),
    "eta": abs(params_tuned[2]),
    "zeta": abs(params_tuned[3])
}
print("Predicted cymatic params (tuned):", cymatic_params)

# -------------------------------
# Generate nodal pattern
# -------------------------------
Nx, Ny = 300, 300
x = np.linspace(0, 1, Nx)
y = np.linspace(0, 1, Ny)
Xg, Yg = np.meshgrid(x, y)

Z = np.sin(cymatic_params["alpha"] * np.pi * Xg) * np.sin(cymatic_params["beta"] * np.pi * Yg) \
    + cymatic_params["eta"] * np.cos(cymatic_params["zeta"] * np.pi * (Xg + Yg))

# Show nodal lines (where Z ~ 0)
plt.figure(figsize=(6,6))
plt.contour(Xg, Yg, Z, levels=[0], colors='black')
plt.title(f"Cymatic Pattern at {f_query} Hz")
plt.axis('off')
plt.show()