From cf0c85947cf4ca814d5d303dc6f5f1bd813967d1 Mon Sep 17 00:00:00 2001
From: weisse02 <andreas.weissenbrunner@ptb.de>
Date: Mon, 27 Feb 2023 12:56:45 +0100
Subject: [PATCH] uncertainty calcul for only one parameterset, cleaned pod
 results

---
 GUI_Elbow.py | 39 ++++++++++++++++++++++++++-------------
 1 file changed, 26 insertions(+), 13 deletions(-)

diff --git a/GUI_Elbow.py b/GUI_Elbow.py
index 01defca..57ce6f8 100644
--- a/GUI_Elbow.py
+++ b/GUI_Elbow.py
@@ -627,6 +627,7 @@ class mainPanel(wx.Panel):
         
         if self.pathint.any():
             if calc_all:
+                print("calculating uncertainty, this will take some hours..")
                 if flow.case == "SingleElbow":
                     self.fm_mean = np.zeros((len(flow.Rk),len(flow.dist),len(flow.phi)))
                     self.fm_std  = np.zeros((len(flow.Rk),len(flow.dist),len(flow.phi)))
@@ -641,8 +642,6 @@ class mainPanel(wx.Panel):
                     self.regint_mean = interpol.RegularGridInterpolator((flow.Rk,flow.dist,flow.phi[:,0]),self.fm_mean,bounds_error = False,fill_value=None)
                     self.regint_std  = interpol.RegularGridInterpolator((flow.Rk,flow.dist,flow.phi[:,0]),self.fm_std,bounds_error = False,fill_value=None)
                 else:
-                    self.fm_mean = np.zeros((len(flow.Rk),len(flow.dl),len(flow.dist),len(flow.phi)))
-                    self.fm_std  = np.zeros((len(flow.Rk),len(flow.dl),len(flow.dist),len(flow.phi)))
                     for i,rc in enumerate(flow.Rk):
                         for j,dl in enumerate(flow.dl):
                             for k,dist in enumerate(flow.dist):
@@ -657,17 +656,24 @@ class mainPanel(wx.Panel):
                     self.regint_std  = interpol.RegularGridInterpolator((flow.Rk,flow.dl,flow.dist,flow.phi[:,0]),self.fm_std,bounds_error = False,fill_value=None)
                 #
             # the uncertainty is only calculated for the one case that is selected
+            # and the other selected parameters 
             else:
                 self.fm_mean = np.zeros((len(flow.dist),len(flow.phi)))
                 self.fm_std  =  np.zeros((len(flow.dist),len(flow.phi)))
-                print("calculating uncertainty for distance:")
-                for j,dist in enumerate(flow.dist):
-                    print(dist)
-                    for k,phi in enumerate(flow.phi[:,0]):
-                        self.fm_mean[j,k], self.fm_std[j,k] = self.get_uncertainty(self.Rc,self.dl,dist,phi)
+                print("calculating uncertainty..")
+                print("for the parameters:")
+                print(str([self.Rc,self.dl,self.dist,self.phi]))
+                print("with the unceraintie:")
+                print(str([self.Rc_u,self.dl_u,self.dist_u,self.phi_u]))
+                fm_mean, fm_std = self.get_uncertainty(self.Rc,self.dl,self.dist,self.phi)
                 #
-                self.regint_mean = interpol.RegularGridInterpolator((flow.dist,flow.phi[:,0]),self.fm_mean,bounds_error = False,fill_value=None)
-                self.regint_std  = interpol.RegularGridInterpolator((flow.dist,flow.phi[:,0]),self.fm_std,bounds_error = False,fill_value=None)
+                print("Expected value u_dis_path:")
+                print(fm_mean)
+                print("Standard uncertainty u_dis_path:")
+                print(fm_std)
+                
+                #self.regint_mean = interpol.RegularGridInterpolator((flow.dist,flow.phi[:,0]),self.fm_mean,bounds_error = False,fill_value=None)
+                #self.regint_std  = interpol.RegularGridInterpolator((flow.dist,flow.phi[:,0]),self.fm_std,bounds_error = False,fill_value=None)
             self.draw_meter()
         #
     def interpol_meter(self,Rc,dl,dist,phi):
@@ -993,9 +999,13 @@ class mainPanel(wx.Panel):
         # np.tensordot(a,ii,axes=0).T
     #
     def get_uncertainty(self,Rc,dl,dist,phi):
-        # Rc = flow.Rk[-1]
-        # test 
+        # performs a quasi-montecarlo for all parameters
+        # this takes a while
+        # it uses interpolation on the flow meter values not on the 
+        # velocities
+        #
         # collect integration values for dist and rc and phi
+        # number of samples in each  dimension
         Ndist = 10 #2* self.bounds_phi[]/(flow.Nphi-1)
         Nrc   = 10
         Nphi  = 81
@@ -1031,7 +1041,7 @@ class mainPanel(wx.Panel):
         #
         if flow.case ==  "SingleElbow":
             # gather all 
-            rcgrid, distgrid, phigrid = np.meshgrid(all_Rc,all_dist,all_phi)
+            rcgrid, distgrid, phigrid = np.meshgrid(all_Rc,all_dist,all_phi,indexing='ij')
             # ugrid = flow.get_profile(distgrid,rcgrid)
             all_fm = self.regint_fm((rcgrid,distgrid,phigrid))
             # flatten the first two dimensions
@@ -1046,9 +1056,12 @@ class mainPanel(wx.Panel):
                 all_dl = dl
             # and ghater all
             
-            rcgrid,dlgrid, distgrid, phigrid = np.meshgrid(all_Rc,all_dl,all_dist,all_phi)
+            rcgrid,dlgrid, distgrid, phigrid = np.meshgrid(all_Rc,all_dl,all_dist,all_phi,indexing='ij')
             # ugrid = flow.get_profile(distgrid,rcgrid)
+            print("discretisation: " + str(rcgrid.shape))
+            print("start interpolation")
             all_fm = self.regint_fm((rcgrid,dlgrid,distgrid,phigrid))
+            print("end interpolation")
         return np.mean(all_fm), np.std(all_fm)
         
     #
-- 
GitLab