import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import rcParams

# Example: set default font globally (optional)
rcParams['font.family'] = 'monospace'      # e.g., 'serif', 'sans-serif', 'monospace'
rcParams['font.size'] = 6

# Load the data: 100 rows, 3 columns
data = np.loadtxt("data.txt")

x, y, z = data[:,0], data[:,1], data[:,2]

limit=1

# Reshape into 10x10 grids
X = x.reshape((101, 101))
Y = y.reshape((101, 101))
Z = z.reshape((101, 101))

def draw(limit,file):
    # Create two masks without modifying Z
    Z_high = np.where(Z >= limit, Z, np.nan)   # terrain cmap
    Z_low  = np.where(Z < limit, Z, np.nan)    # blue

    Z_sea_level = np.where(Z >= 0, Z, np.nan) 
    
    global_count = Z.size
    count = np.sum(Z < limit)
    pourcentage = (count / global_count) * 100
    
    print(f"{limit}m  {count} {pourcentage:5.2f}")
    
    # 3D surface plot
    fig = plt.figure(figsize=(8, 6), dpi=300)
    ax = fig.add_subplot(111, projection="3d")

    surf1 = ax.plot_surface(X, Y, Z_high, cmap="terrain")
    surf2 = ax.plot_surface(X, Y, Z_low, color="blue")

    # Optional: top view
    # ax.view_init(elev=50, azim=-90)
    # ax.view_init(elev=90, azim=-90)
    ax.view_init(elev=90, azim=-90)
    
    ax.set_title("Elevation "+str(limit)+"m", fontsize=16, fontweight="bold")
    #ax.set_xlabel("X")
    #ax.set_ylabel("Y")
    #ax.set_zlabel("Z")
    ax.set_xticks([])
    ax.set_yticks([])
    ax.set_zticks([])
    #plt.show()
    plt.savefig(file)

for i in range(1,15):
    draw(i,"alt"+str(i)+".png")
