import os
import matplotlib.pyplot as plt
import matplotlib.patches as patches
import numpy as np

# Stations data (in meters)
stations_x = np.array([0, 1167, 2333, 3500, 4667, 5833, 7000]) / 1000.0
z_keel = np.array([306.3, 109.5, 17.8, 0.0, 50.9, 279.7, 1341.7]) / 1000.0
z_sheer = np.array([1201.7, 1129.1, 1065.6, 1050.0, 1099.3, 1214.6, 1341.7]) / 1000.0

mb_25 = np.array([707.9, 908.9, 1012.5, 1050.2, 912.8, 515.5, 0.0]) / 1000.0
mb_50 = np.array([827.7, 1056.5, 1164.8, 1201.8, 1070.8, 645.0, 0.0]) / 1000.0
mb_75 = np.array([868.5, 1105.3, 1212.5, 1248.0, 1125.3, 702.0, 0.0]) / 1000.0
mb_deck = np.array([875.0, 1112.7, 1219.2, 1254.2, 1134.3, 716.2, 0.0]) / 1000.0

# Interpolator function for sections using piecewise cubic / linear
def get_station_section(st_idx):
    zk = z_keel[st_idx]
    zs = z_sheer[st_idx]
    if st_idx == 6:
        return lambda z: np.zeros_like(z)
    z_lev = np.array([zk, zk + 0.25*(zs-zk), zk + 0.50*(zs-zk), zk + 0.75*(zs-zk), zs])
    y_lev = np.array([0.0, mb_25[st_idx], mb_50[st_idx], mb_75[st_idx], mb_deck[st_idx]])
    
    def interp_fn(z_query):
        z_q = np.asarray(z_query)
        res = np.interp(z_q, z_lev, y_lev, left=0.0, right=y_lev[-1])
        # Smooth rounding at bilge
        return res
    return interp_fn

section_interps = [get_station_section(i) for i in range(7)]

# Waterlines (WLs)
wls_z = [0.150, 0.250, 0.350, 0.600, 0.850, 1.050]
wl_y_at_stations = {}
for wz in wls_z:
    y_vals = []
    for i in range(7):
        zk = z_keel[i]
        zs = z_sheer[i]
        if wz < zk or wz > zs or i == 6:
            y_vals.append(0.0)
        else:
            y_vals.append(float(section_interps[i](wz)))
    wl_y_at_stations[wz] = np.array(y_vals)

# Setup Drawing Canvas (Estilo Brana / Técnico Naval)
plt.style.use('default')
fig = plt.figure(figsize=(19, 13), facecolor='#FFFFFF')
gs = fig.add_gridspec(3, 2, height_ratios=[1.15, 1.15, 1.25], width_ratios=[2.3, 1.0], hspace=0.32, wspace=0.20)

ax_profile = fig.add_subplot(gs[0, 0])
ax_plan = fig.add_subplot(gs[1, 0])
ax_body = fig.add_subplot(gs[2, 0])
ax_info = fig.add_subplot(gs[:, 1])

# =========================================================================
# 1. PLANO DE PERFIL (SHEER PLAN) COM ESPELHO DE POPA A 20° E CALADO 350 MM
# =========================================================================
ax_profile.set_facecolor('#FCFDFE')
ax_profile.set_title("1. PLANO DE PERFIL (SHEER PLAN) — ESPELHO DE POPA INCLINADO A 20° & DWL = 350 mm", 
                     fontsize=11, fontweight='bold', color='#1A365D', pad=8)

# Linha d'Água de Projeto (DWL = 350 mm)
dwl_z = 0.350
ax_profile.axhline(dwl_z, color='#2B6CB0', linestyle='--', linewidth=1.4, label="Linha d'Água de Projeto (DWL = 350 mm)")
ax_profile.fill_between([-0.7, 7.5], dwl_z, -0.15, color='#EBF8FF', alpha=0.55)

# Waterlines lines
for wz in [0.150, 0.250, 0.600, 0.850, 1.050]:
    ax_profile.axhline(wz, color='#CBD5E0', linestyle=':', linewidth=0.7)
    ax_profile.text(-0.55, wz, f"WL {int(wz*1000)}", color='#718096', fontsize=6.5, va='center')

ax_profile.text(-0.55, dwl_z, "DWL 350", color='#2B6CB0', fontsize=7.5, fontweight='bold', va='center')

# Espelho de popa inclinado a 20° para ré
# Base na Baliza 1 (X=0.0, Z=0.3063 m), Topo em X = -0.3259 m, Z = 1.2017 m
z_transom = np.linspace(z_keel[0], z_sheer[0], 50)
angle_rad = np.radians(20.0)
x_transom = -(z_transom - z_keel[0]) * np.tan(angle_rad)

# Polynomial curves along hull
poly_sheer = np.poly1d(np.polyfit(stations_x, z_sheer, 3))
poly_keel  = np.poly1d(np.polyfit(stations_x, z_keel, 4))

x_hull = np.linspace(0.0, 7.0, 200)
z_sheer_fine = poly_sheer(x_hull)
z_keel_fine  = np.clip(poly_keel(x_hull), 0, None)
z_keel_fine[-1] = z_sheer_fine[-1]

# Plot profile lines
ax_profile.plot(x_transom, z_transom, color='#1A365D', linewidth=2.4, label="Espelho de Popa (20° de inclinação)")
ax_profile.plot(x_hull, z_sheer_fine, color='#1A365D', linewidth=2.4, label="Linha de Borda / Tosamento")
ax_profile.plot(x_hull, z_keel_fine, color='#1A365D', linewidth=2.4, label="Linha de Quilha / Roda de Proa")

# Shading hull area
ax_profile.fill_between(x_hull, z_keel_fine, z_sheer_fine, color='#CBD5E0', alpha=0.35)

# Transom angle indicator arc
arc = patches.Arc((0.0, z_keel[0]), 0.6, 0.6, angle=90, theta1=0, theta2=20, color='#C53030', linewidth=1.5)
ax_profile.add_patch(arc)
ax_profile.text(-0.16, z_keel[0]+0.22, "20°", color='#C53030', fontweight='bold', fontsize=8.5)
ax_profile.plot([0.0, 0.0], [z_keel[0], z_sheer[0]], color='#A0AEC0', linestyle='-.', linewidth=0.9)

# Balizas verticais
for i, st in enumerate(stations_x):
    ax_profile.axvline(st, color='#A0AEC0', linestyle=':', linewidth=0.8)
    ax_profile.text(st, -0.08, f"Bal. {i+1}\n({int(st*1000)} mm)", color='#4A5568', fontsize=7.2, ha='center')

# Cotas LOA
ax_profile.annotate("", xy=(-0.326, 1.55), xytext=(7.0, 1.55), arrowprops=dict(arrowstyle="<->", color="#1A365D", lw=1.2))
ax_profile.text(3.34, 1.63, "LOA = 7,33 m (Com inclinação do espelho) | LPP = 7,00 m (Balizas 1 a 7)", 
                color='#1A365D', fontweight='bold', fontsize=8.5, ha='center')

# Calado Cota
ax_profile.annotate("", xy=(3.5, 0.0), xytext=(3.5, dwl_z), arrowprops=dict(arrowstyle="<->", color="#2B6CB0", lw=1.3))
ax_profile.text(3.65, 0.175, "Calado T = 350 mm", color='#2B6CB0', fontweight='bold', fontsize=8, va='center')

ax_profile.set_xlim(-0.7, 7.6)
ax_profile.set_ylim(-0.18, 1.85)
ax_profile.set_ylabel("Altura Acima da Base Z (m)", fontsize=8.5, fontweight='bold', color='#1A365D')
ax_profile.legend(loc='upper right', fontsize=7.5, framealpha=0.92)
ax_profile.grid(True, linestyle=':', alpha=0.45)

# =========================================================================
# 2. PLANTA DE MEIAS-BOCAS E LINHAS D'ÁGUA (HALF-BREADTH PLAN)
# =========================================================================
ax_plan.set_facecolor('#FCFDFE')
ax_plan.set_title("2. PLANO DE LINHAS D'ÁGUA (HALF-BREADTH PLAN) — MEIAS-BOCAS & CONVÉS", 
                  fontsize=11, fontweight='bold', color='#1A365D', pad=8)
ax_plan.axhline(0, color='#E53E3E', linestyle='-.', linewidth=1.3, label="Linha de Centro (CL)")

poly_deck = np.poly1d(np.polyfit(stations_x, mb_deck, 4))
y_deck_fine = np.clip(poly_deck(x_hull), 0, None)
y_deck_fine[-1] = 0.0

ax_plan.plot(x_hull, y_deck_fine, color='#1A365D', linewidth=2.4, label="Linha de Borda (Boca Máx = 2.508 mm)")
ax_plan.plot(x_hull, -y_deck_fine, color='#1A365D', linewidth=2.4)

# DWL 350 mm
poly_dwl = np.poly1d(np.polyfit(stations_x, wl_y_at_stations[0.350], 4))
y_dwl_fine = np.clip(poly_dwl(x_hull), 0, None)
y_dwl_fine[-1] = 0.0

ax_plan.plot(x_hull, y_dwl_fine, color='#2B6CB0', linestyle='--', linewidth=1.6, label="DWL 350 mm (Linha d'Água de Projeto)")
ax_plan.plot(x_hull, -y_dwl_fine, color='#2B6CB0', linestyle='--', linewidth=1.6)

# Other WLs
colors_wl = {0.150: '#4FD1C5', 0.250: '#319795', 0.600: '#805AD5', 0.850: '#DD6B20', 1.050: '#718096'}
for wz in [0.150, 0.250, 0.600, 0.850]:
    poly_wz = np.poly1d(np.polyfit(stations_x, wl_y_at_stations[wz], 4))
    y_wz_fine = np.clip(poly_wz(x_hull), 0, None)
    y_wz_fine[-1] = 0.0
    ax_plan.plot(x_hull, y_wz_fine, color=colors_wl[wz], linestyle=':', linewidth=1.1, label=f"WL {int(wz*1000)} mm")
    ax_plan.plot(x_hull, -y_wz_fine, color=colors_wl[wz], linestyle=':', linewidth=1.1)

# Espelho de popa fechamento
ax_plan.plot([0, 0], [-mb_deck[0], mb_deck[0]], color='#1A365D', linewidth=2.4)
ax_plan.fill_between(x_hull, -y_deck_fine, y_deck_fine, color='#E2E8F0', alpha=0.35)

# Balizas
for i, (st, y_val) in enumerate(zip(stations_x, mb_deck)):
    ax_plan.plot([st, st], [-y_val, y_val], color='#A0AEC0', linestyle=':', linewidth=0.8)
    ax_plan.plot(st, y_val, 'o', color='#2B6CB0', markersize=3.5)
    ax_plan.plot(st, -y_val, 'o', color='#2B6CB0', markersize=3.5)
    ax_plan.text(st, -1.52, f"Bal. {i+1}", color='#4A5568', fontsize=7.2, ha='center')

# Cota Boca Máxima
ax_plan.annotate("", xy=(3.5, -1.254), xytext=(3.5, 1.254), arrowprops=dict(arrowstyle="<->", color="#C53030", lw=1.4))
ax_plan.text(3.65, 0.0, "Boca Máx = 2.508 mm (Baliza 4)", color='#C53030', fontweight='bold', fontsize=8, va='center')

ax_plan.set_xlim(-0.7, 7.6)
ax_plan.set_ylim(-1.75, 1.75)
ax_plan.set_xlabel("Comprimento Longitudinal X (m) [Balizas 1 a 7]", fontsize=8.5, fontweight='bold', color='#1A365D')
ax_plan.set_ylabel("Semi-Boca Y (m)", fontsize=8.5, fontweight='bold', color='#1A365D')
ax_plan.legend(loc='upper right', fontsize=7.2, framealpha=0.92, ncol=2)
ax_plan.grid(True, linestyle=':', alpha=0.45)

# =========================================================================
# 3. PLANO DE BALIZAS (BODY PLAN) — SEÇÕES EXATAS DA TABELA
# =========================================================================
ax_body.set_facecolor('#FCFDFE')
ax_body.set_title("3. PLANO DE BALIZAS (BODY PLAN) — SEÇÕES TRANSVERSAIS EXATAS DA TABELA", 
                  fontsize=11, fontweight='bold', color='#1A365D', pad=8)
ax_body.axvline(0, color='#E53E3E', linestyle='-.', linewidth=1.3)
ax_body.axhline(dwl_z, color='#2B6CB0', linestyle='--', linewidth=1.3, label="DWL 350 mm")

colors_st = ['#2D3748', '#4A5568', '#718096', '#E53E3E', '#3182CE', '#00B4D8', '#38A169']

for i in range(7):
    zk = z_keel[i]
    zs = z_sheer[i]
    if i == 6:
        continue
    z_pts = np.linspace(zk, zs, 60)
    y_pts = section_interps[i](z_pts)
    
    if i < 3: # Balizas 1, 2, 3 (Popa / Ré - Lado Esquerdo)
        ax_body.plot(-y_pts, z_pts, color=colors_st[i], linewidth=2.0 if i==0 else 1.5, label=f"Bal. {i+1} ({int(stations_x[i]*1000)} mm)")
        ax_body.text(-y_pts[-1]-0.04, z_pts[-1], f"Bal {i+1}", fontsize=7.5, fontweight='bold', color=colors_st[i], ha='right')
        z_dots = np.array([zk, zk + 0.25*(zs-zk), zk + 0.50*(zs-zk), zk + 0.75*(zs-zk), zs])
        y_dots = np.array([0.0, mb_25[i], mb_50[i], mb_75[i], mb_deck[i]])
        ax_body.plot(-y_dots, z_dots, 'o', color=colors_st[i], markersize=3)
    else: # Balizas 4, 5, 6 (Proa / Vante - Lado Direito)
        ax_body.plot(y_pts, z_pts, color=colors_st[i], linewidth=2.2 if i==3 else 1.5, label=f"Bal. {i+1} ({int(stations_x[i]*1000)} mm)")
        ax_body.text(y_pts[-1]+0.04, z_pts[-1], f"Bal {i+1}", fontsize=7.5, fontweight='bold', color=colors_st[i], ha='left')
        z_dots = np.array([zk, zk + 0.25*(zs-zk), zk + 0.50*(zs-zk), zk + 0.75*(zs-zk), zs])
        y_dots = np.array([0.0, mb_25[i], mb_50[i], mb_75[i], mb_deck[i]])
        ax_body.plot(y_dots, z_dots, 'o', color=colors_st[i], markersize=3)

# Waterlines reference lines on Body Plan
for wz in wls_z:
    ax_body.axhline(wz, color='#CBD5E0', linestyle=':', linewidth=0.7)

ax_body.text(-0.75, 0.08, "SEÇÕES DE POPA / RÉ\n(Balizas 1, 2, 3)", color='#4A5568', fontsize=8, fontweight='bold', ha='center')
ax_body.text(0.75, 0.08, "SEÇÕES DE PROA / VANTE\n(Balizas 4, 5, 6)", color='#3182CE', fontsize=8, fontweight='bold', ha='center')

ax_body.set_xlim(-1.55, 1.55)
ax_body.set_ylim(-0.12, 1.55)
ax_body.set_xlabel("Semi-Boca Y (m) [Bombordo (-) / Boreste (+)]", fontsize=8.5, fontweight='bold', color='#1A365D')
ax_body.set_ylabel("Altura Z (m)", fontsize=8.5, fontweight='bold', color='#1A365D')
ax_body.legend(loc='lower left', fontsize=7.0, framealpha=0.92, ncol=2)
ax_body.grid(True, linestyle=':', alpha=0.45)

# =========================================================================
# 4. PAINEL DE DADOS OFICIAIS & TABELA DE COTAS BRANA
# =========================================================================
ax_info.axis('off')
ax_info.set_facecolor('#FCFDFE')

# Hydrostatic calculation at T = 350 mm
station_areas_dwl = []
for i in range(7):
    zk = z_keel[i]
    zs = z_sheer[i]
    if dwl_z <= zk or i == 6:
        station_areas_dwl.append(0.0)
    else:
        z_sub = np.linspace(zk, dwl_z, 40)
        y_sub = section_interps[i](z_sub)
        area_half = np.trapezoid(y_sub, z_sub)
        station_areas_dwl.append(2.0 * area_half)

# Simpson integration for volume
dx_st = 1.1667
vol_dwl = (dx_st / 3.0) * (station_areas_dwl[0] + 4*station_areas_dwl[1] + 2*station_areas_dwl[2] + 
                           4*station_areas_dwl[3] + 2*station_areas_dwl[4] + 4*station_areas_dwl[5] + station_areas_dwl[6])
disp_dwl_ton = vol_dwl * 1.025 # Água salgada BTS
cb_dwl = vol_dwl / (7.0 * 2.508 * 0.350)
am_dwl = station_areas_dwl[3]
cm_dwl = am_dwl / (2.508 * 0.350)
cp_dwl = vol_dwl / (7.0 * am_dwl)

info_text = (
    "MEMORIAL DO PLANO DE LINHAS OFICIAL\n"
    "=================================================\n"
    "Embarcação: Saveiro Ciclone 2 (7,0 m)\n"
    "Projetista / Aluno: Rafael Rosa Batista\n"
    "Orientação: Profa. Magdalena (Mada) - Brana / EEE\n"
    "Grupo: Grupo 8 | Malha Base: 291,67 mm\n"
    "-------------------------------------------------\n"
    "PARÂMETROS GEOMÉTRICOS PRINCIPAIS:\n"
    "• Comprimento Entre Balizas (LPP): 7,000 m\n"
    "• Comprimento Total (LOA): 7,326 m (c/ espelho 20°)\n"
    "• Boca Moldada Máxima: 2,508 m (Baliza 4)\n"
    "• Pontal Moldado (D): 1,050 m (Meia-Nau)\n"
    "• Calado de Projeto (T): 0,350 m (350 mm)\n"
    "• Borda Livre de Projeto (F): 0,700 m (700 mm)\n"
    "• Inclinação do Espelho de Popa: 20,0° para ré\n"
    "-------------------------------------------------\n"
    "HIDROSTÁTICA OFICIAL (DWL = 350 mm / Água Salgada):\n"
    f"• Volume de Carena (V): {vol_dwl:.3f} m³\n"
    f"• Deslocamento em Carga (Δ): {disp_dwl_ton*1000:.1f} kg ({disp_dwl_ton:.2f} t)\n"
    f"• Área da Seção Mestra (Am): {am_dwl:.3f} m²\n"
    f"• Coeficiente de Bloco (Cb): {cb_dwl:.3f}\n"
    f"• Coeficiente de Seção Mestra (Cm): {cm_dwl:.3f}\n"
    f"• Coeficiente Prismático (Cp): {cp_dwl:.3f}\n"
    "-------------------------------------------------\n"
    "TABELA DE COTAS (mm):\n"
    "Stn   X      Z_Keel  Z_Borda  MB_25   MB_50   MB_Deck\n"
    " 1    0.0    306.3   1201.7   707.9   827.7    875.0\n"
    " 2  1167.0   109.5   1129.1   908.9  1056.5   1112.7\n"
    " 3  2333.0    17.8   1065.6  1012.5  1164.8   1219.2\n"
    " 4  3500.0     0.0   1050.0  1050.2  1201.8   1254.2\n"
    " 5  4667.0    50.9   1099.3   912.8  1070.8   1134.3\n"
    " 6  5833.0   279.7   1214.6   515.5   645.0    716.2\n"
    " 7  7000.0  1341.7   1341.7     0.0     0.0      0.0\n"
    "================================================="
)

ax_info.text(0.02, 0.98, info_text, transform=ax_info.transAxes,
             fontsize=7.3, fontfamily='monospace', verticalalignment='top',
             bbox=dict(boxstyle='round,pad=0.55', facecolor='#FFFFFF', edgecolor='#CBD5E0', linewidth=1.2))

out_img = "/home/magdalena/eee/alunos/rafael-rosa-batista/Plano_de_Linhas_Oficial_Ciclone_2_Calado_350mm.png"
plt.savefig(out_img, dpi=300, bbox_inches='tight')
plt.close()
print("High-res technical plan generated successfully at:", out_img)