Added plot bars for DBP, MBP, SBP, Q

functions.py

@@ -5,6 +5,7 @@ Created on Thu Jul 16 15:08:27 2020
 @author: Aloma Blanch
 """
 
+import matplotlib
 import numpy as np
 import pandas as pd
 import matplotlib.pyplot as plt
@@ -47,7 +48,7 @@ def cycle(folder,dt,save_path):
     ax.legend(loc=0)
     # plt.show()
     plt.savefig(save_path + '/Pressure_max_n_cycles.pdf')
-    return (floor(mean(peak_tdiff)*1000)/1000,len(t_pk_loc))
+    return (floor(mean(peak_tdiff[1:])*1000)/1000,len(t_pk_loc))
 
 def error_plot(folder,t_step,r_criteria,save_path):
     # Load iterations and residual error
@@ -105,7 +106,7 @@ def periodicity(project,folder,dt,T_cyc,n_cyc,save_path):
     
     fig, ax = plt.subplots()
     ax.plot(time,pressure[:,-1]/1333.22,'b')
-    ax.plot(time[peak_P_pos], peak_P,'ro',label='Cylce pike')
+    ax.plot(time[peak_P_pos], peak_P,'ro',label='Cycle pike')
     ax.set(xlabel='time [s]', ylabel='Pressure [mmHg]',
         title='Pressure at last outlet')
     ax.spines['right'].set_visible(False)
@@ -219,4 +220,68 @@ def rcr(project_folder):
         Rc_C_Rd[i][2] = float(rcr_lines[4+i*6][:-1])
     return Rc_C_Rd
         
-        
+def barPlot(project_folder,DBP,MBP,SBP,PP,Q_avg,save_path):
+    import csv
+    with open(project_folder+'/aimed.csv', 'r') as f:
+        aimed_all = list(csv.reader(f, delimiter=';'))  
+        
+    def plot_bar(CFD,aimed,name,save_path,unit):
+        labels = []
+        for i in range(0,len(CFD)):
+            labels.append('ROI-'+str(i+2))
+        x = np.arange(len(labels))  # the label locations
+        width = 0.38  # the width of the bars
+        fig, ax = plt.subplots()
+        rects1 = ax.bar(x - width/2, CFD, width, label='CFD',color='k')
+        rects2 = ax.bar(x + width/2, aimed, width, label='Aimed',color='white', edgecolor='black', hatch='/')
+        ax.set_ylabel(name + unit)
+        ax.set_title(name)
+        ax.set_xticks(x)
+        ax.set_xticklabels(labels)
+        ax.legend()
+        high = max(max(aimed,CFD))
+        ax.set_ylim(top = 1.7*high)
+        def autolabel(rects):
+            """Attach a text label above each bar in *rects*, displaying its height."""
+            for rect in rects:
+                height = rect.get_height()
+                ax.annotate('{}'.format(height),
+                            xy=(rect.get_x() + rect.get_width() / 2, height),
+                            xytext=(0,3),  # 3 points vertical offset
+                            textcoords="offset points",
+                            ha='center', va='bottom',
+                            fontsize = 7.5)
+    
+        autolabel(rects1)
+        autolabel(rects2)
+        fig.tight_layout()
+        # plt.show()
+        plt.savefig(save_path + '/'+name+'_bar.pdf')
+        fig.savefig(save_path + '/'+name+'_bar.jpg',dpi=150)        
+    
+    # DBP
+    CFD = [round(float(i),1) for i in DBP]
+    aimed = [round(float(i),1) for i in aimed_all[0]]
+    plot_bar(CFD,aimed,'DBP',save_path,' [mmHg]')
+    
+    # MBP
+    CFD = [round(float(i),1) for i in MBP]
+    aimed = [round(float(i),1) for i in aimed_all[1]]
+    plot_bar(CFD,aimed,'MBP',save_path,' [mmHg]')
+    
+    # SBP
+    CFD = [round(float(i),1) for i in SBP]
+    aimed = [round(float(i),1) for i in aimed_all[2]]
+    plot_bar(CFD,aimed,'SBP',save_path,' [mmHg]')
+    
+    # PP
+    CFD = [round(float(i),1) for i in PP]
+    aimed = [round(float(i),1) for i in aimed_all[3]]
+    plot_bar(CFD,aimed,'PP',save_path,' [mmHg]')
+    
+    # DBP
+    CFD = [round(float(i),1) for i in Q_avg]
+    aimed = [round(float(i),1) for i in aimed_all[4]]
+    plot_bar(CFD,aimed,'Q_avg',save_path,' [mL/s]')
+        
+    

generatePDF.py

@@ -135,14 +135,23 @@ def generatePDF(save_path,project,DBP,MBP,SBP,PP,Q_avg,txt1,txt2,Rc_C_Rd,T_cyc,n
         pdf.cell(width_cell[0],10,'ROI-'+ str(i+2),1,0,'C')
         pdf.cell(width_cell[1],10,str(round(Q_avg[i],3)),1,1,'C')
 
+    pdf.print_section(4, 'Checking CFD values vs aimed')
+    pdf.image(save_path +'/DBP_bar.jpg', x = 40, y = 50, w = 140, h = 100, type = '', link = '')
+    pdf.image(save_path +'/MBP_bar.jpg', x = 40, y = 170, w = 140, h = 100, type = '', link = '')
+    pdf.print_section(4, 'Checking CFD values vs aimed')
+    pdf.image(save_path +'/SBP_bar.jpg', x = 40, y = 50, w = 140, h = 100, type = '', link = '')
+    pdf.image(save_path +'/PP_bar.jpg', x = 40, y = 170, w = 140, h = 100, type = '', link = '')
+    pdf.print_section(4, 'Checking CFD values vs aimed')
+    pdf.image(save_path +'/Q_avg_bar.jpg', x = 40, y = 50, w = 140, h = 100, type = '', link = '')
+    pdf.cell(0, 10,'', 0, 1)
     
-    pdf.print_section(4, 'Checking Inlet Flow Waveform and time steps saved')
+    pdf.print_section(5, 'Checking Inlet Flow Waveform and time steps saved')
     pdf.image(save_path +'/inlet_waveform.jpg', x = 40, y = 50, w = 140, h = 100, type = '', link = '')
     pdf.set_font('Arial', '', 10)
     pdf.cell(0, 5,txt2, 0, 1)
     
     width_cell=[30,20,110];
-    pdf.print_section(5, 'Outlet Boundary Condition Applied')
+    pdf.print_section(6, 'Outlet Boundary Condition Applied')
     pdf.set_font('Arial', '', 10)
     pdf.cell(0, 10,'3-Element Windkessel model outlet boundary condition applied in SimVascular', 0, 1)
     pdf.cell(0, 5,' ', 0, 1)
@@ -162,7 +171,7 @@ def generatePDF(save_path,project,DBP,MBP,SBP,PP,Q_avg,txt1,txt2,Rc_C_Rd,T_cyc,n
         pdf.cell(width_cell[0],10,str(Rc_C_Rd[i][1]),1,0,'C') # Third column of row 1 
         pdf.cell(width_cell[0],10,str(Rc_C_Rd[i][2]),1,1,'C') # Third column of row 1
 
-    pdf.print_section(6, 'Solver Parameters')
+    pdf.print_section(7, 'Solver Parameters')
     pdf.set_font('Arial', '', 10)
     pdf.cell(0, 10,'Cardiac cycle period - ' + str(T_cyc) + 's', 0, 1)
     pdf.cell(0, 10,'Number of cycles - ' + str(n_cyc), 0, 1)

main.py

@@ -20,7 +20,7 @@ import os.path
 # from fpdf import FPDF
 
 
-from functions import error_plot, periodicity, pressure, flow, inlet_flow_waveform, rcr, cycle
+from functions import error_plot, periodicity, pressure, flow, inlet_flow_waveform, rcr, cycle, barPlot
 from generatePDF import generatePDF  
 
 
@@ -53,7 +53,7 @@ Rc_C_Rd = rcr(project_folder)
 
 # Extracting number of cycles and period of cardiac cycle
 (T_cyc,n_cyc) = cycle(folder,dt,save_path)
-        
+T_cyc=0.477    
 # Cehcking convergency and periodicity
 error_plot(folder,dt,rc,save_path)
 txt1 = periodicity(project,folder,dt,T_cyc,n_cyc,save_path)
@@ -67,5 +67,8 @@ txt1 = periodicity(project,folder,dt,T_cyc,n_cyc,save_path)
 # Inlet Flow Rate + and t saved
 txt2 = inlet_flow_waveform(project_folder,t_btw_rst,N_ts,dt,T_cyc,n_cyc,save_path)
 
+# Extract bar plots CFD vs. aimed
+barPlot(project_folder,DBP,MBP,SBP,PP,Q_avg,save_path)
+
 # Create PDF report
 generatePDF(save_path,project,DBP,MBP,SBP,PP,Q_avg,txt1,txt2,Rc_C_Rd,T_cyc,n_cyc,N_ts,dt)