From 0a949ffd3d8ff37c9387ff236eeff67be6ae06a3 Mon Sep 17 00:00:00 2001
From: Aloma Blanch <aloma_b@hotmail.com>
Date: Thu, 16 Jul 2020 20:02:02 -0500
Subject: [PATCH] adding pressure plot

---
 functions.py | 19 ++++++++++++++++++-
 main.py      |  4 +++-
 2 files changed, 21 insertions(+), 2 deletions(-)

diff --git a/functions.py b/functions.py
index 886d970..562dc55 100644
--- a/functions.py
+++ b/functions.py
@@ -42,7 +42,7 @@ def error_plot(folder,t_step,r_criteria,save):
     # Semilog scale
     plt.figure()  
     plt.semilogy(time,error)
-    plt.semilogy(time,r_criteria*np.ones(len(error)),'r:')
+    plt.semilogy(time,r_criteria*np.ones(len(error)),'r')
     plt.ylabel('Residual error')
     plt.xlabel('Time steps')
     plt.title('Log - Last nonlinear residual error for each time step')
@@ -76,3 +76,20 @@ def periodicity(project,folder,dt,T_cyc,n_cyc):
     if (peak_Pdiff[-1]<=1): print('The numerical simulation \'{0}\' has achieve periodicity!\nSystolic Blood Pressure (SBP):\nsecond-last cycle = {1:.2f} mmHg,\nlast cycle = {2:.2f} mmHg,\n\u0394mmHg = {3:.2f} mmHg'.format(project,peak_P[-2],peak_P[-1],peak_Pdiff[-1]))  
 
 
+def pressure(folder,N_ts,T_cyc,dt,n_cyc):
+    pressure = np.loadtxt(folder+'/PHistRCR.dat',skiprows=2,)
+    Nc = round(T_cyc/dt)
+    time = np.linspace(0,T_cyc,Nc)
+    fig, ax = plt.subplots()
+    for i in range(0,pressure.shape[1]):
+        ax.plot(time,pressure[N_ts-Nc:N_ts,-i]/1333.22,label='ROI-'+str(i+2))
+        
+    ax.set(xlabel='time [s]', ylabel='Pressure [mmHg]',
+        title='Pressure @ each outlet')
+    ax.spines['right'].set_visible(False)
+    ax.spines['top'].set_visible(False)
+    ax.legend(loc=0)
+    plt.show()
+        
+        
+        
\ No newline at end of file
diff --git a/main.py b/main.py
index d48278c..ead6b2a 100644
--- a/main.py
+++ b/main.py
@@ -13,7 +13,7 @@ from scipy import signal
 import statistics
 
 
-from functions import error_plot, periodicity
+from functions import error_plot, periodicity, pressure
   
 # Selct dir 
 Tk().withdraw()
@@ -48,4 +48,6 @@ rc = float(mylines[26][18:-1])
 error_plot(folder,dt,rc,False)
 periodicity(project,folder,dt,T_cyc,n_cyc)
 
+# Pressure
+pressure(folder,N_ts,T_cyc,dt,n_cyc)