quickplot_amountgas.py

import os
import re
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

path = "C:\\Users\\User\\Documents\\Documents\\University\\Bachelor_RUG\\Year_4\\IIB\\Research_project\\Research\\Simulations\\"
files = os.listdir(path)
files.remove("Analysis")
lst_volatiles1 = ["H2O", "CO2", "H2", "CH4", "CO", "N2", "NH3", "O2"]
lst_volatiles2 = ["H2O", "CO2", "H2", "CH4", "CO", "N2", "NH3"]

for i in files:
    print(i)
    module = i.split("_")[1]
    if module == "JANUS":
        lst_volatiles = lst_volatiles1
    else:
        lst_volatiles = lst_volatiles2

    fig, ax1 = plt.subplots(1,1)
    x_gas = dict.fromkeys(lst_volatiles)
    p = dict.fromkeys(lst_volatiles)
    pp_gas = dict.fromkeys(lst_volatiles)
    for w in x_gas.keys():
        x_gas[w] = []
        p[w] = []
        pp_gas[w] = []


    for t in lst_volatiles:
        data = os.listdir(path+i+f"\\Results\\Data\\{t}\\")

        A = []
        for k in data:
            a = int(re.search(r"\d+", k).group())
            A.append(a)
        A = sorted(A)                          #And pull time points from the names of the files
        Time = A

        A_new = []
        for k in A:
            a = "phasestte"+str(k)+".txt"
            A_new.append(a)

        data = A_new

        for n in data:
            df = pd.read_csv(path+i+f"\\Results\\Data\\{t}\\{n}", header = None, sep=" ")
            df.columns = df.iloc[0]
            df = df[1:]
            l = len(df["x_gas"])-2
            x_gas[t].append(float(df["x_gas"][l])) #sum([float(w) for w in df["x_gas"]]) / len([float(w) for w in df["x_gas"]]) ) 
            p[t].append(float(df["Pressure(Pa)"][l])) #sum([float(w) for w in df["Pressure(Pa)"]]) / len([float(w) for w in df["Pressure(Pa)"]]) )
        
        x_gas[t] = np.array(x_gas[t])
        p[t] = np.array(p[t])/1e5
        pp_gas[t] = x_gas[t]*p[t]

        ax1.plot(Time, pp_gas[t], linewidth = 1, label = f"{t}")
        ax1.set_xscale("log")
        ax1.set_yscale("log")
        ax1.set_xlabel("Time (yr)")
        ax1.set_ylabel("Partial pressure of volatiles at the surface (bar)")
        ax1.legend(loc='lower left', bbox_to_anchor=(1, 0), borderaxespad=0.2)

    fig.savefig(f"C:\\Users\\User\\Documents\\Documents\\University\\Bachelor_RUG\\Year_4\\IIB\\Research_project\\Research\\Simulations\\{i}\\Amountgas.png", bbox_inches = "tight", dpi=150)