TFTB_GUIproject/Sun_process/Magnetic_Field.py

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from datetime import date, timedelta
def find_max_min(arr):
    #print(arr)
    if max(arr) > abs(min(arr)):
        return arr / max(arr), max(arr / max(arr))

    else:
        return arr / abs(min(arr)), abs(min(arr / abs(min(arr))))

def elipse(kx, ky, kz):
    u = np.linspace(0, 2 * np.pi, 20)
    v = np.linspace(0, np.pi, 20)

    return kx * np.outer(np.cos(u), np.sin(v)), ky * np.outer(np.sin(u), np.sin(v)), kz * np.outer(np.ones(np.size(u)), np.cos(v))

def Magnetic_3d_picture(x, y, z, begin, end, step, df_B, count, numb):

    for i in range(1, step + 1):
        # x1, kx = find_max_min(x[begin + (i - 1) * count:begin + i * count])
        # y1, ky = find_max_min(y[begin + (i - 1) * count:begin + i * count])
        # z1, kz = find_max_min(z[begin + (i - 1) * count:begin + i * count])

        x1, kx = find_max_min(x[begin + (i - 1) * count:begin + i * count])
        y1, ky = find_max_min(y[begin + (i - 1) * count:begin + i * count])
        z1, kz = find_max_min(z[begin + (i - 1) * count:begin + i * count])

        x2, y2, z2 = elipse(kx, ky, kz)
        with open ("datatime.csv", "a") as file:
            pd.set_option('display.max_rows', None)
            file.write(f"{x1} {y1} {z1}\n")
        fig = plt.figure(figsize=(10, 10))
        ax = fig.add_subplot(projection='3d')
        ax.scatter(x1, y1, z1, c='black', marker=">", s = 10)
        ax.plot_wireframe(x2, y2, z2, linewidth=0.1, color='b')
        date_begin = date(2023, 1, 1) + timedelta(days=df_B.iloc[begin + (i - 1) * count]["Day"]-1)
        date_end = date(2023, 1, 1) + timedelta(days=df_B.iloc[begin + (i) * count]["Day"] - 1)
        plt.title("Days: " + str(date_begin) +" : " + str(date_end) + " Hours: " + str(df_B.iloc[begin + (i - 1) * count]["Hour"]) + " - " + str(df_B.iloc[begin - 1 + (i) * count]["Hour"]), size=20,  weight="bold")

        #plt.suptitle(str(count)+"  Days №"+str(df_B.iloc[begin + (i - 1) * count]["Day"]) + " - " + str(df_B.iloc[begin + (i) * count]["Day"] ))
        plt.xlabel('$Bx$', size=32, color='black', labelpad=25)
        ax.tick_params(axis="x", labelsize=24)
        plt.ylabel('$By$', size=32, color='black', labelpad=25)
        ax.tick_params(axis="y", labelsize=24)
        ax.set_zlabel('$Bz$', size=32, labelpad=25)
        ax.tick_params(axis="z", labelsize=24)
        #ax.view_init(30, -126) # изменение угла 3d графика
        #plt.savefig(f"/home/kmg/Pictures/Graphs/diagram_{numb}.png", dpi=300)
        plt.show()

def Read_B_Table(file):
    #arr = ["Year", "Day", "Hour", "Minute", "Sec", "Bx_GSE", "By_GSE", "Bz_GSE"]
    arr = ["Year", "Day", "Hour", "Minute", "Bx_GSE", "By_GSE", "Bz_GSE"]
    df = pd.read_table(file, names=arr, sep="\s+")
    for i in arr:
       df.loc[(df[i] > 9999)] = 0
    #print(df.shape[0])
    return df

def SMA(df1, df2, df3, window):
    rolling_mean_1 = pd.Series(df1).rolling(window=window).mean()
    rolling_mean_2 = pd.Series(df2).rolling(window=window).mean()
    rolling_mean_3 = pd.Series(df3).rolling(window=window).mean()

    fig = plt.figure(figsize=(12,12), dpi=100)
    ax1=fig.add_subplot(311)
    ax1.plot(df1, label='AMD 50 Day SMA', color='black')
    ax1.plot(rolling_mean_1, label='AMD 50 Day SMA', color='orange')
    ax1.plot([0 for i in range(len(df1))], color = 'magenta')
    plt.xlabel("Bx")

    ax2 = fig.add_subplot(312)
    ax2.plot(df2, label='AMD 50 Day SMA', color='black')
    ax2.plot(rolling_mean_2, label='AMD 50 Day SMA', color='orange')
    ax2.plot([0 for i in range(len(df2))], color='magenta')
    plt.xlabel("By")

    ax3 = fig.add_subplot(313)
    ax3.plot(df3, label='AMD 50 Day SMA', color='black')
    ax3.plot(rolling_mean_3, label='AMD 50 Day SMA', color='orange')
    ax3.plot([0 for i in range(len(df3))], color='magenta')
    plt.xlabel("Bz")
    plt.show()
    nuli = []
    for i, j in enumerate(rolling_mean_1): # находит переходы через ноль
        if rolling_mean_1[i] < 0 and rolling_mean_1[i - 1] > 0:
            print("из - в +",i, j)
            nuli.append(i)
        elif rolling_mean_1[i] > 0 and rolling_mean_1[i-1] < 0:
            print("из + в  -",i, j)
            nuli.append(i)

    return rolling_mean_1, rolling_mean_2, rolling_mean_3

def plot_SMA (smaX,smaY, smaZ, begin, end):

    fig = plt.figure(figsize=(10, 10))
    ax = fig.add_subplot(projection='3d')
    arr1 = []
    arr2 = []

    for i in range(begin, end):
        if smaX[i] > 0:
            arr1.append(i)
        else:
            arr2.append(i)

    ax.scatter(smaX[arr1], smaY[arr1], smaZ[arr1], c='black', marker=">", s=10)
    ax.scatter(smaX[arr2], smaY[arr2], smaZ[arr2], c='red', marker="<", s=10)

    plt.xlabel('$Bx$', size=24, color='black', labelpad=20)
    plt.ylabel('$By$', size=24, color='black', labelpad=20)
    ax.set_zlabel('$Bz$', size=24, labelpad=20)
    plt.show()

    return arr1, arr2

def Magnetic_3d_picture_sma(x, y, z, begin, end, step, df_B, count, numb, arr1, arr2):

    for i in range(1, step + 1):
        # x1, kx = find_max_min(x[begin + (i - 1) * count:begin + i * count])
        # y1, ky = find_max_min(y[begin + (i - 1) * count:begin + i * count])
        # z1, kz = find_max_min(z[begin + (i - 1) * count:begin + i * count])

        x1, kx = find_max_min(x[begin + (i - 1) * count:begin + i * count])
        y1, ky = find_max_min(y[begin + (i - 1) * count:begin + i * count])
        z1, kz = find_max_min(z[begin + (i - 1) * count:begin + i * count])

        x2, y2, z2 = elipse(kx, ky, kz)
        with open ("datatime.csv", "a") as file:
            pd.set_option('display.max_rows', None)
            file.write(f"{x1} {y1} {z1}\n")
        fig = plt.figure(figsize=(10, 10))
        ax = fig.add_subplot(projection='3d')

        ax.scatter(x1[arr1], y1[arr1], z1[arr1], c='black', marker=">", s=10)
        ax.scatter(x1[arr2], y1[arr2], z1[arr2], c='red', marker=">", s=10)

        ax.plot_wireframe(x2, y2, z2, linewidth=0.1, color='b')
        date_begin = date(2023, 1, 1) + timedelta(days=df_B.iloc[begin + (i - 1) * count]["Day"]-1)
        date_end = date(2023, 1, 1) + timedelta(days=df_B.iloc[begin + (i) * count]["Day"] - 1)
        plt.title("Days: " + str(date_begin) +" : " + str(date_end) + " Hours: " + str(df_B.iloc[begin + (i - 1) * count]["Hour"]) + " - " + str(df_B.iloc[begin - 1 + (i) * count]["Hour"]), size=20,  weight="bold")

        #plt.suptitle(str(count)+"  Days №"+str(df_B.iloc[begin + (i - 1) * count]["Day"]) + " - " + str(df_B.iloc[begin + (i) * count]["Day"] ))
        plt.xlabel('$Bx$', size=32, color='black', labelpad=25)
        ax.tick_params(axis="x", labelsize=24)
        plt.ylabel('$By$', size=32, color='black', labelpad=25)
        ax.tick_params(axis="y", labelsize=24)
        ax.set_zlabel('$Bz$', size=32, labelpad=25)
        ax.tick_params(axis="z", labelsize=24)
        #ax.view_init(30, -126) # изменение угла 3d графика
        #plt.savefig(f"/home/kmg/Pictures/Graphs/diagram_{numb}.png", dpi=300)
        plt.show()

def B_Diagrams(df, begin, end, numb):
    X = df["Bx_GSE"]
    Y = df["By_GSE"]
    Z = df["Bz_GSE"]

    count = end - begin
    #count = 5000
    # #720 # шаг по 3 часа для 15сек
    step = (end - begin) // count

    window = 7000 # окно для скользящей средней
    smaX, smaY, smaZ=SMA(X, Y, Z, window)
    arr1, arr2 = plot_SMA(smaX, smaY, smaZ, begin, end)

    Magnetic_3d_picture_sma(X, Y, Z, begin, end, step, df, count, numb, arr1, arr2)

    Magnetic_3d_picture(X, Y, Z, begin, end, step, df, count, numb)
    #Magnetic_3d_picture(smaX, smaY, smaZ, begin, end, step, df, count, numb)
    # combination_rotation_animation(X[begin:end], Y[begin:end], Z[begin:end], begin, step, 720)
    # Anime(X, Y, Z, step)

if __name__ == "__main__":
    df_test = Read_B_Table("sun.lst")
    print()