From 50e9a7300f5b66bbdeaf77af47c6cb9c8fd15884 Mon Sep 17 00:00:00 2001
From: farid <farid.mohammadi@iws.uni-stuttgart.de>
Date: Mon, 17 May 2021 17:18:58 +0200
Subject: [PATCH] [PA-A] modified post plot.

---
 .../tests/PA-A/ffpm_validation_stokesdarcy.py |  4 +-
 .../tests/PA-A/ffpm_validation_stokespnm.py   |  2 +-
 BayesValidRox/tests/PA-A/post_plot.py         | 54 +++++++++++--------
 3 files changed, 35 insertions(+), 25 deletions(-)

diff --git a/BayesValidRox/tests/PA-A/ffpm_validation_stokesdarcy.py b/BayesValidRox/tests/PA-A/ffpm_validation_stokesdarcy.py
index 1ff7de7a4..4cba2ebe9 100644
--- a/BayesValidRox/tests/PA-A/ffpm_validation_stokesdarcy.py
+++ b/BayesValidRox/tests/PA-A/ffpm_validation_stokesdarcy.py
@@ -130,7 +130,7 @@ def run(params, errorPerc=0.1, couplingcond='BJ', PCEEDMethod='normal'):
     MetaModelOpts = Metamodel(Inputs)
     
     # Select if you want to preserve the spatial/temporal depencencies
-    MetaModelOpts.DimRedMethod = 'PCA'
+    # MetaModelOpts.DimRedMethod = 'PCA'
     # MetaModelOpts.varPCAThreshold = 99.99
     
     # Select your metamodel method
@@ -296,7 +296,7 @@ def run(params, errorPerc=0.1, couplingcond='BJ', PCEEDMethod='normal'):
     # Select the inference method
     BayesOptsCalib.SamplingMethod = 'MCMC'
     # BayesOptsCalib.MultiProcessMCMC = False
-    BayesOptsCalib.MCMCnwalkers = 300 #100
+    BayesOptsCalib.MCMCnwalkers = 500 #100
     # BayesOptsCalib.MCMCnSteps = 1500 #5000
     
     BayesOptsCalib.MAP ='mean'
diff --git a/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py b/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py
index b06dd8ec3..db8222f18 100755
--- a/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py
+++ b/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py
@@ -291,7 +291,7 @@ def run(params, errorPerc=0.1, PCEEDMethod='normal'):
     BayesOptsCalib.SamplingMethod = 'MCMC'
     # import emcee
     # BayesOptsCalib.MCMCmoves = emcee.moves.StretchMove()
-    BayesOptsCalib.MCMCnwalkers = 300 #300
+    BayesOptsCalib.MCMCnwalkers = 500 #300
     # BayesOptsCalib.MCMCnSteps = 1500 #5000
     
     BayesOptsCalib.MAP ='mean'
diff --git a/BayesValidRox/tests/PA-A/post_plot.py b/BayesValidRox/tests/PA-A/post_plot.py
index de4d77a9f..9f26ebec2 100755
--- a/BayesValidRox/tests/PA-A/post_plot.py
+++ b/BayesValidRox/tests/PA-A/post_plot.py
@@ -6,7 +6,7 @@ Created on Tue Feb 18 08:28:13 2020
 @author: farid
 """
 import numpy as np
-import os
+import os, sys
 import seaborn as sns
 import h5py
 from scipy.stats import norm
@@ -20,7 +20,8 @@ from matplotlib.patches import Patch
 import matplotlib.lines as mlines
 from matplotlib import ticker
 import matplotlib.pylab as plt
-
+# Add BayesValidRox path
+# sys.path.insert(0,'./../../../BayesValidRox/')
 plt.style.use('seaborn-deep') # talk, paper, poster
 fsize = 50
 params = {'legend.fontsize': fsize,
@@ -34,7 +35,7 @@ params = {'legend.fontsize': fsize,
 #          'mathtext.fontset' : 'stix',
 #          'mathtext.rm'      : 'serif',
           'font.family'      : 'serif',
-          'font.serif'       : "Times New Roman", # or "Times" 
+          'font.serif'       : "Times New Roman", # or "Times"
           'figure.figsize'   : (32, 24),
           # 'savefig.dpi':1500,
           'text.usetex': True
@@ -50,22 +51,30 @@ def upper_rugplot(data, height=.05, ax=None, **kwargs):
                     axis=-1)
     lc = LineCollection(segs, transform=ax.get_xaxis_transform(), **kwargs)
     ax.add_collection(lc)
-    
+
 def postPredictiveplot(modelName, errorPrec, case='Calib'):
     # result_folder = './Results_10_03_2021/outputs_{}/'.format(modelName.split('-v')[0])
     result_folder = './'
     directory = result_folder+'Outputs_Bayes_'+modelName+'_'+case
     OutputDir = ('postPred_'+modelName+'_'+case)
     if not os.path.exists(OutputDir): os.makedirs(OutputDir)
-    
+
     data = pd.read_csv('data/stokesData'+case+'.csv')
-    
-    
-    
+
+    # Load Post pred file
     f = h5py.File(directory+'/'+"postPredictive.hdf5", 'r+')
+
+    # Load PCEModel
+    with open(result_folder+'PCEModel_'+modelName+'.pkl', 'rb') as input:
+        PCEModel = pickle.load(input)
+
+    # Generate Prior Predictive
+    #priorPred, std = PCEModel.eval_metamodel(nsamples=10000)
+
+
     for OutputName in ['velocity [m/s]', 'p']:
         # x_coords = np.array(f["x_values/"+OutputName])
-        
+
         # Find pointIDs
         csv_file = 'pressure_points.csv' if OutputName == 'p' else 'velocity_points.csv'
         if case == 'Calib':
@@ -77,29 +86,30 @@ def postPredictiveplot(modelName, errorPrec, case='Calib'):
 
         for idx, x in enumerate(pointIDs):
             fig, ax = plt.subplots(1, 1)
-            
-            
+
+            # Prior predictive
+            #outputs = priorPred[OutputName][:,idx]
+            #sns.histplot(outputs[outputs>0], ax=ax,# bins=50,
+                         #color='blue', alpha=0.4,stat="count")
+
             # Posterior predictive
             postPred = np.array(f["EDY/"+OutputName])[:,idx]
             sns.histplot(postPred[postPred>0], ax=ax,# bins=50,
                          color='orange',stat="count") #normalizes counts so that the sum of the bar heights is 1
-            
+
             # Reference data from the pore-scale simulation
             # ax.axvline(x=data[OutputName][idx], linewidth=15, color='green')
             sns.histplot(np.random.normal(data[OutputName][idx], errorPrec*data[OutputName][idx], len(postPred[postPred>0])), ax=ax,
                           color='green', alpha=0.5, stat="count")
-            
+
             # Print conidence interval of ExpDesign.Y (Trained area)
-            with open(result_folder+'PCEModel_'+modelName+'.pkl', 'rb') as input:
-                PCEModel = pickle.load(input)
-            
             modelRuns = PCEModel.ExpDesign.Y[OutputName][:,idx]
             upper_rugplot(modelRuns, ax=ax, alpha=0.75, color='grey')
-            
+
             # fmt = ticker.StrMethodFormatter("{x}")
             # ax.xaxis.set_major_formatter(fmt)
             # ax.yaxis.set_major_formatter(fmt)
-            
+
             legend_elements = [
                 Patch(facecolor='orange', edgecolor='orange',label='Posterior Pred.'),
                 Patch(facecolor='green', edgecolor='green',alpha=0.5, label='Ref. Data'),
@@ -107,7 +117,7 @@ def postPredictiveplot(modelName, errorPrec, case='Calib'):
                               linestyle='None', markersize=75, markeredgewidth=1.5, label='Orig. Responses')]
 
             ax.legend(handles=legend_elements, fontsize=75)
-            
+
             font = {'family': 'serif',
                     'weight': 'normal',
                     'size': 100,
@@ -131,6 +141,6 @@ def postPredictiveplot(modelName, errorPrec, case='Calib'):
             fig.savefig('./'+OutputDir+'/PointID_'+str(idx+1)+'_'+plotname+'.pdf', bbox_inches='tight')
             plt.close()
 
-# modelName = 'ffpm-stokespnm' #stokespnm stokesdarcyER stokesdarcyBJ
-# postPredictiveplot(modelName, case='Calib')
-# postPredictiveplot(modelName+'-valid', case='Valid')
+# modelName = 'ffpm-stokesdarcyER' #stokespnm stokesdarcyER stokesdarcyBJ
+# postPredictiveplot(modelName, errorPrec=0.1, case='Calib')
+# postPredictiveplot(modelName+'-valid', errorPrec=0.1, case='Valid')
-- 
GitLab