import matplotlib.pyplot as plt
import numpy as np
import integral_calculation as ic


def vel_field_plot(vinf, alpha, x_pan, y_pan, phi, panel_len, lambda_list):

    # Grid Parameters
    x_pts = 30
    y_pts = 30
    x_range = [-0.3, 1.5]
    y_range = [-0.4, 0.6]

    # Generate Grid
    x_grid = np.linspace(x_range[0], x_range[1], x_pts)
    y_grid = np.linspace(y_range[0], y_range[1], y_pts)
    X, Y = np.meshgrid(x_grid, y_grid)

    # Initialize velocities
    u_vel = np.zeros([len(x_grid), len(y_grid)])
    v_vel = np.zeros([len(x_grid), len(y_grid)])
    tot_vel = np.zeros([len(x_grid), len(y_grid)])

    # Compute panel's number
    panel_num = len(x_pan) - 1

    # For loop to compute velocities in each point
    for x in range(len(x_grid)):
        for y in range(len(y_grid)):

            # Compute integrals for 1 point
            # Mx, My = ic.compute_Mpj(x_pan, y_pan, X[x, y], Y[x, y], phi, panel_len)
            Mx, My = teste.STREAMLINE_SPM(X[x,y],Y[x,y],x_pan,y_pan,phi,panel_len)
            mx_sum = 0
            my_sum = 0

            # Sum integral's values for each panel
            for j in range(panel_num):
                mx_sum += (lambda_list[j]/(2*np.pi)) * Mx[j]
                my_sum += (lambda_list[j]/(2*np.pi)) * My[j]

            # Compute velocities
            u_vel[x,y] = (vinf * np.cos(alpha)) + mx_sum
            v_vel[x,y] = (vinf * np.sin(alpha)) + my_sum
            tot_vel[x,y] = (u_vel[x,y]**2 + v_vel[x,y]**2) ** 0.5

    ##%% Plotting

    # Plot airfoil
    plt.fill(x_pan, y_pan, 'grey')
    plt.axis('equal')

    # Plot total velocity vector field
    contf = plt.contourf(X, Y, tot_vel, levels=np.linspace(0, 1.6*vinf ,100), extend='both')
    cbar = plt.colorbar(contf,orientation='horizontal',shrink=0.5, pad = 0.1)
    cbar.set_label('Velocity', fontsize=16)
    plt.title('Contour of velocity field', fontsize=16);
    plt.show()


def vel_streamlines(vinf, alpha, x_pan, y_pan, phi, panel_len, lambda_list):

    # Grid Parameters
    x_pts = 40
    y_pts = 40
    x_range = [-0.3, 1.5]
    y_range = [-0.4, 0.6]

    # Generate Grid
    x_grid = np.linspace(x_range[0], x_range[1], x_pts)
    y_grid = np.linspace(y_range[0], y_range[1], y_pts)
    X, Y = np.meshgrid(x_grid, y_grid)

    # Initialize velocities
    u_vel = np.zeros([len(x_grid), len(y_grid)])
    v_vel = np.zeros([len(x_grid), len(y_grid)])

    # Compute panel's number
    panel_num = len(x_pan) - 1

    # For loop to compute velocities in each point
    for x in range(len(x_grid)):
        for y in range(len(y_grid)):

            # Compute integrals for 1 point
            Mx, My = ic.compute_Mpj(x_pan, y_pan, X[x, y], Y[x, y], phi, panel_len)
            mx_sum = 0
            my_sum = 0

            # Sum integral's values for each panel
            for j in range(panel_num):
                mx_sum += (lambda_list[j]/(2*np.pi)) * Mx[j]
                my_sum += (lambda_list[j]/(2*np.pi)) * My[j]

            # Compute velocities
            u_vel[x,y] = (vinf * np.cos(alpha)) + mx_sum
            v_vel[x,y] = (vinf * np.sin(alpha)) + my_sum

    print(u_vel)
    # Plot airfoil
    plt.fill(x_pan, y_pan, 'grey')
    plt.axis('equal')
    plt.quiver(X, Y,u_vel, v_vel)
    plt.show()


def panel_vis(x_pan, y_pan, xc, yc, panel_len, phi):

    # Compute panel's angle with x axis
    delta = phi + np.pi/2

   # Create list with normal points
    panel_num = len(x_pan) - 1
    norm_pts_1 = list()
    norm_pts_2 = list()

   # Compute points to plot panel's normals
    for j in range(panel_num):
        x1 = xc[j]
        y1 = yc[j]

        x2 = xc[j] + panel_len[j]*np.cos(delta[j])
        y2 = yc[j] + panel_len[j]*np.sin(delta[j])

        norm_pts_1.append([x1, x2])
        norm_pts_2.append([y1, y2])

    ##%% PLOTTING

    #Plot airfoil
    plt.plot(x_pan, y_pan, 'black', linewidth=0.75)

    #Plot Control points
    plt.scatter(list(xc), list(yc))

    #Plot normals
    for n in range(panel_num):
        x_val = norm_pts_1[n]
        y_val = norm_pts_2[n]
        plt.plot(x_val, y_val, 'grey', linewidth=0.75)

    plt.axis('equal')
    plt.show()