118 lines
3.5 KiB
Python
118 lines
3.5 KiB
Python
#!/usr/bin/env python
|
|
import numpy as np
|
|
import matplotlib.pyplot as plt
|
|
|
|
T_cyc = 0.477
|
|
N_ts = 14280
|
|
cyc = 5
|
|
|
|
pressure_filename = "Z:\\Aloma\\Class\\final_BIEN6931\\Simulations\\MCS_m04\\96-procs_case\\PHistRCR.dat"
|
|
#pressure_filename = "Z:\\Aloma\\Class\\final_BIEN6931\\Simulations\\sim03_RCR_2cyc\\PHistRCR_ok.dat"
|
|
#pressure = np.loadtxt(pressure_filename)
|
|
pressure = np.loadtxt(pressure_filename, skiprows=1,)
|
|
|
|
time = np.linspace(0,T_cyc*cyc,pressure.shape[0])
|
|
|
|
# Ploting Pressure over iterations
|
|
# use this to have the plot in different window %matplotlib qt
|
|
n=2
|
|
plt.figure()
|
|
plt.plot(time,pressure[:,n]/1333.22)
|
|
plt.ylabel('Pressure [mmHg]')
|
|
plt.xlabel('Time [s]')
|
|
plt.title('Presure evolution ROI-2')
|
|
axes = plt.gca()
|
|
xmax = time[-1]
|
|
xmin = time[0]
|
|
axes.set_xlim([xmin,xmax])
|
|
ymin = min(pressure[:,n]/1333.22)
|
|
ymax = max(pressure[:,n]/1333.22+1)
|
|
axes.set_ylim([ymin,ymax])
|
|
plt.show()
|
|
|
|
# Ploting Pressure in all the vessels same plot over time
|
|
n = pressure.shape
|
|
i=0
|
|
plt.figure()
|
|
while (i<n[1]):
|
|
plt.plot(time[round(N_ts-N_ts/cyc):N_ts],pressure[round(N_ts-N_ts/cyc):N_ts,i]/1333.22, label='ROI-'+str(i+2))
|
|
i=i+1
|
|
|
|
plt.ylabel('Pressure [mmHg]')
|
|
plt.xlabel('Time [s]')
|
|
#plt.title('5th Cardiac Cylce')
|
|
axes = plt.gca()
|
|
#xmax = time[-1]
|
|
#xmin = time[round(N_ts-N_ts/cyc)]
|
|
#axes.set_xlim([xmin,xmax])
|
|
#ymin = np.amin(pressure[round(N_ts-N_ts/cyc):N_ts,:]/1333.22)
|
|
#ymax = np.amax(pressure[round(N_ts-N_ts/cyc):N_ts,:]/1333.22)+1
|
|
#axes.set_ylim([ymin,ymax])
|
|
plt.legend(loc='upper right')
|
|
plt.show()
|
|
|
|
|
|
# Ploting Pressure over time
|
|
i=0
|
|
while (i<n[1]):
|
|
plt.plot(time,pressure[:,i]/1333.22)
|
|
#plt.plot(round(N_ts/cyc), flowRate[round(N_ts/cyc),i], 'ro')
|
|
plt.ylabel('Pressure [mmHg]')
|
|
plt.xlabel('Time [s]')
|
|
plt.title('Pressure waveform ROI-'+str(i+2)+' [cgs]')
|
|
#axes = plt.gca()
|
|
xmax = time[-1]
|
|
xmin = time[0]
|
|
axes.set_xlim([xmin,xmax])
|
|
ymin = min(pressure[:,0]/1333.22)
|
|
ymax = max(pressure[:,0]/1333.22)
|
|
axes.set_ylim([ymin,ymax])
|
|
plt.show()
|
|
i = i + 1
|
|
|
|
|
|
# Plot Flow Rate to analyze that inlet is as expected all outlets all cycles
|
|
flowRate_filename = "Z:\\Aloma\\Class\\final_BIEN6931\\Simulations\\MCS_m04\\96-procs_case\\QHistRCR.dat"
|
|
#flowRate_filename = "Z:\\Aloma\\Class\\final_BIEN6931\\Simulations\\sim03_RCR_cyc\\QHistRCR_ok.dat"
|
|
#flowRate = np.loadtxt(flowRate_filename)
|
|
flowRate = np.loadtxt(flowRate_filename, skiprows=1,)
|
|
|
|
n = flowRate.shape
|
|
i=0
|
|
|
|
while (i<n[1]):
|
|
plt.plot(time,flowRate[:,i])
|
|
#plt.plot(round(N_ts/cyc), flowRate[round(N_ts/cyc),i], 'ro')
|
|
plt.ylabel('Flow Rate Q [mL/s]')
|
|
plt.xlabel('Time [s]')
|
|
plt.title('Outlet Flow Rate ROI-'+str(i+2)+ ' [cgs]')
|
|
#axes = plt.gca()
|
|
xmax = time[-1]
|
|
xmin = time[0]
|
|
axes.set_xlim([xmin,xmax])
|
|
ymin = min(flowRate[:,0])
|
|
ymax = max(flowRate[:,0])
|
|
axes.set_ylim([ymin,ymax])
|
|
plt.show()
|
|
i = i + 1
|
|
|
|
|
|
# Modifications need to be done in order to Plot only the cycle that user's want
|
|
cycle = int(input("Which cycle you whant to plot (velocity plot)? "))
|
|
n = flowRate.shape
|
|
i=0
|
|
while (i<n[1]):
|
|
plt.plot(time[round(N_ts-N_ts/cyc):N_ts],flowRate[round(N_ts-N_ts/cyc):N_ts,i])
|
|
i = i + 1
|
|
|
|
plt.ylabel('Flow Rate Q [mL/s]')
|
|
plt.xlabel('Time [s]')
|
|
plt.title('Out Flow Rate cylce '+str(cycle)+ '- all ROI [cgs]')
|
|
#axes = plt.gca()
|
|
xmax = time[-1]
|
|
xmin = time[0]
|
|
axes.set_xlim([xmin,xmax])
|
|
ymin = min(flowRate[:,0])
|
|
ymax = max(flowRate[:,0])
|
|
axes.set_ylim([ymin,ymax])
|
|
plt.show() |