Changeset 25836 for issm/trunk/src/m/qmu/dakota_out_parse.py

Timestamp:

12/08/20 08:45:53 (4 years ago)

Author:

Mathieu Morlighem

Message:

merged trunk-jpl and trunk for revision 25834

Location:

issm/trunk

Files:

• . (modified) (1 prop)
• src (modified) (1 prop)
• src/m/qmu/dakota_out_parse.py (modified) (26 diffs)

issm/trunk/src/m/qmu/dakota_out_parse.py

@@ -8 +8 @@
 from helpers import *
 
-#Note: this may be re-written later to take advantage of Python's file i / o mechanics
-#    as it is written now it is often difficult to work with, but is analagous to
-#    the Matlab version of dakota_out_parse
+# NOTE: May be rewritten later to take advantage of Python's file I/O 
+# mechanics. As it is written now, it is often difficult to work with, but is 
+# analagous to the MATLAB version of dakota_out_parse.
 
 
 def dakota_out_parse(filei):  # {{{
-    '''
-  read a Dakota .out or .dat output file and parse it.
-
-  [method, dresp, scm, pcm, srcm, prcm] = dakota_out_parse(filei)
-
-  where the required input is:
-    filei         (character, name of .out file)
-
-  the required output is:
-    method        (character, dakota method name)
-    dresp         (structure array, responses)
-
-  and the optional output is:
-    scm           (double array, simple correlation matrix)
-    pcm           (double array, partial correlation matrix)
-    srcm          (double array, simple rank correlation matrix)
-    prcm          (double array, partial rank correlation matrix)
-
-  the filei will be prompted if empty.  the fields of dresp
-  are particular to the data contained within the file.  the
-  scm, pcm, srcm, and prcm are output by dakota only for the
-  sampling methods.
-
-  this function reads a dakota .out output file and parses it
-  into the matlab workspace.  it operates in a content - driven
-  fashion, where it skips the intermediate data and then parses
-  whatever output data it encounters in the order in which it
-  exists in the file, rather than searching for data based on
-  the particular method.  (this makes it independent of method.)
-  it also can read and parse the .dat tabular_output file.
-
-  this data would typically be used for plotting and other
-  post - processing within matlab or excel.
-'''
+    """DAKOTA_OUT_PARSE - read a Dakota .out or .dat output file and parse it.
+
+    Usage:
+        [method, dresp, scm, pcm, srcm, prcm] = dakota_out_parse(filei)
+
+    where the required input is,
+        filei         (character, name of .out file)
+
+    the required output is,
+        method        (character, Dakota method name)
+        dresp         (structure array, responses)
+
+    and the optional output is,
+        scm           (double array, simple correlation matrix)
+        pcm           (double array, partial correlation matrix)
+        srcm          (double array, simple rank correlation matrix)
+        prcm          (double array, partial rank correlation matrix)
+
+    The filei will be prompted for if empty. The fields of dresp are particular 
+    to the data contained within the file. The scm, pcm, srcm, and prcm are 
+    output by Dakota only for the sampling methods.
+
+    This function reads a Dakota .out output file and parses it into the Python 
+    runtime. It operates in a content-driven fashion, where it skips the 
+    intermediate data and then parses whatever output data it encounters in the 
+    order in which it exists in the file, rather than searching for data based 
+    on the particular method (this makes it independent of method). It also can 
+    read and parse the .dat tabular_output file.
+
+    This data would typically be used for plotting and other postprocessing 
+    within MATLAB or Excel.
+
+    TODO:
+    - Figure out why output from Dakota is different under MATLAB and Python 
+    (is it the input file that we write?)
+
+    "Copyright 2009, by the California Institute of Technology. ALL RIGHTS 
+    RESERVED. United States Government Sponsorship acknowledged. Any commercial 
+    use must be negotiated with the Office of Technology Transfer at the 
+    California Institute of Technology. (NTR 47078)
+
+    This software may be subject to U.S. export control laws. By accepting this 
+    software, the user agrees to comply with all applicable U.S. export laws 
+    and regulations. User has the responsibility to obtain export licenses, or 
+    other export authority as may be required before exporting such information 
+    to foreign countries or providing access to foreign persons."
+    """
+
     if filei is None:
         help(dakota_out_parse)
@@ -53 +67 @@
 
     if not isfile(filei) or getsize(filei) == 0:
-        filei = str(eval(input('Input file?  ')))
+        filei = str(eval(input('Input file? ')))
 
     #fidi = fopen(sprintf('%s', filei), 'r')
@@ -62 +76 @@
         fline = fidi.readline()
         if getsize(filei) == 0 or fline == '':
-            raise RuntimeError('File ' + filei + ' is empty.')
-
-        dresp = []  # of struct()
+            raise RuntimeError('File ' + filei + ' is empty')
+
+        dresp = [] # of struct()
         scm = struct()
         pcm = struct()
@@ -77 +91 @@
             fidi.seek(0, 0)
 
-    #  loop through the file to find the Dakota method name
+        # loop through the file to find the Dakota method name
         fline = findline(fidi, 'method', True)
         if fline is None:
@@ -87 +101 @@
                 [ntokens, tokens] = fltokens(fline)
                 method = tokens[0].strip()
-                print('Dakota method =\'' + method + '\'.')
+                print('Dakota method = \'' + method + '\'')
             elif fline[6] in ['N', 'n']:
                 fline = findline(fidi, 'methodName = ')
                 [ntokens, tokens] = fltokens(fline)
                 method = tokens[2].strip()
-                print('Dakota methodName = "' + method + '".')
-
-    #  loop through the file to find the function evaluation summary
+                print('Dakota methodName = \'' + method + '\'')
+
+        # loop through the file to find the function evaluation summary
         counter = 0
         fline = ''
         nfeval = nfeval_read(fidi, fline)
 
-        #  process each results section based on content of the file
+        # process each results section based on content of the file
         while counter < 10:
-            # because python makes file i / o difficult
+            # because python makes file I/O difficult
             # if we see 10 + blank lines in a row then we have reached EOF
             # (tests show actual maximum number of blank lines is around 5)
@@ -108 +122 @@
             else:
                 counter = 0
-    #     ipos = ftell(fidi)
+            # ipos = ftell(fidi)
             fline = fidi.readline()
             if fline == '' or fline.isspace():
@@ -147 +161 @@
                 'Unexpected line: ' + str(fline)
 
-    #     fidi.seek(ipos, 0)
-
-    #  loop through the file to verify the end
+            # fidi.seek(ipos, 0)
+
+    # loop through the file to verify the end
 
     # fline = findline(fidi, '<<<<< Single Method Strategy completed')
@@ -155 +169 @@
     #     return
     #
-        print('End of file successfully reached.')
+        print('End of file successfully reached')
     #close(fidi)
     #except Exception as err:
     #print "ERROR in dakota_out_parse: " + err
     #raise err
-    #raise RuntimeError(filei + ' could not be opened.')
+    #raise RuntimeError(filei + ' could not be opened')
 
     return [method, dresp, scm, pcm, srcm, prcm]
@@ -167 +181 @@
 
 def dak_tab_out(fidi, fline):  # {{{
-    #  function to parse the dakota tabular output file
-
-    print('Reading Dakota tabular output file.')
-
-    #  process column headings of matrix (skipping eval_id)
+    """DAK_TAB_OUT - function to parse the Dakota tabular output file
+    """
+
+    print('Reading Dakota tabular output file')
+
+    # Process column headings of matrix (skipping eval_id)
     [ntokens, tokens] = fltokens(fline)
 
-    # New file DAKOTA versions > 6
-    if strncmpi(fline, '%eval_id interface', 18):
+    if strncmpi(fline, '%eval_id interface', 18): # Dakota versions >= 6
         offset = 2
-    else:  #DAKOTA versions < 6
+    else:  # Dakota versions < 6
         offset = 1
 
-    desc = [['' for i in range(ntokens - offset)]]
+    desc = ['' for i in range(ntokens - offset)]
     data = np.zeros((1, ntokens - offset))
 
     for i in range(ntokens - offset):
-        desc[0][i] = str(tokens[i + offset])
-
-    print("Number of columns (Dakota V + R)=" + str(ntokens - 2) + '.')
-
-    #  process rows of matrix
+        desc[i] = str(tokens[i + offset])
+
+    print('Number of columns (Dakota V + R) = {}'.format(ntokens - 2))
+
+    # Process rows of matrix
     nrow = 0
     while True:
-
         fline = fidi.readline()
 
@@ -202 +215 @@
         [ntokens, tokens] = fltokens(fline)
 
-    #  add row values to matrix (skipping eval_id)
-
+        # Add row values to matrix (skipping eval_id)
         for i in range(ntokens - offset):
             data[nrow, i] = tokens[i + offset]
@@ -209 +221 @@
         nrow = nrow + 1
 
-    print('Number of rows (Dakota func evals) = ' + str(nrow) + '.')
-
-    #  calculate statistics
-
-    #  since normfit doesn't have a dim argument, and matlab isvector is True
-    #  for a 1xn matrix, handle the case of one row explicitly
-
-    #  Update: normfit_issm.py does handle this case
+    print('Number of rows (Dakota func evals) = ' + str(nrow))
+
+    # Calculate statistics
     if (np.size(data, 0) > 1):
-        #dmean  =mean   (data)
-        #dstddev = std    (data, 0)
+        #dmean = mean(data)
+        #dstddev = std(data, 0)
         [dmean, dstddev, dmeanci, dstddevci] = normfit_issm(data, 0.05)
     else:
@@ -229 +236 @@
             [dmean[0, i], dstddev[0, i], dmeanci[:, i], dstddevci[:, i]] = normfit_issm(data[:, i], 0.05)
 
-    dmin = data.min(0)
+    dmin = data.min(axis=0)
     dquart1 = prctile_issm(data, 25, 0)
-    dmedian = np.median(data, 0)
+    dmedian = np.median(data, axis=0)
     dquart3 = prctile_issm(data, 75, 0)
-    dmax = data.max(0)
+    dmax = data.max(axis=0)
     dmin95 = prctile_issm(data, 5, 0)
     dmax95 = prctile_issm(data, 95, 0)
 
-    # Note: the following line may cause the following warning
-    #    (should not crash or invalidate results) when one of
-    #    the inputs does not change with respect to the
-    #    other / s causing an internal divide-by - zero error
-
-    # / usr / local / lib / python2.7 / dist - packages / numpy / lib / function_base.py:3163:
-    #    RuntimeWarning: invalid value encountered in true_divide
-    #    c / = stddev[:, None]
-
-    #    (and / or the same but with "c / = stddev[None, :]")
-
+    # NOTE: The following line may cause the following warning (should not 
+    # crash or invalidate results) when one of the inputs does not change with 
+    # respect to the other(s), causing an internal divide-by-zero error,
+    #
+    #       /usr/local/lib/python2.7/dist-packages/numpy/lib/function_base.py:3163:
+    #       RuntimeWarning: invalid value encountered in true_divide
+    #       c /= stddev[:, None]
+    #
+    # (and/or the same but with "c /= stddev[None, :]")
+
+    # Equivalent to Dakota scm, MATLAB corrcoef, and Excel correl
     dcorrel = np.corrcoef(data.T)
 
-    #  divide the data into structures for consistency
+    # Divide the data into structures for consistency
     dresp = []
-
     for i in range(len(desc)):
         dresp.append(struct())
-        dresp[i].descriptor = str(desc[i])
+        dresp[i].descriptor = desc[i]
         dresp[i].sample = data[:, i]
         dresp[i].mean = dmean[i]
@@ -294 +300 @@
     [ntokens, tokens] = fltokens(fline)
     nfeval = tokens[4]
-    print('  Dakota function evaluations = ' + str(int(nfeval)) + '.')
+    print('  Dakota function evaluations = ' + str(int(nfeval)))
 
     return nfeval
@@ -309 +315 @@
     [ntokens, tokens] = fltokens(fline)
     nsamp = tokens[3]
-    print('  Dakota samples = ' + str(int(nsamp)) + '.')
+    print('  Dakota samples = ' + str(int(nsamp)))
 
     return nsamp
@@ -340 +346 @@
         dresp[-1].coefvar = tokens[12]
 
-    print('  Number of Dakota response functions = ' + str(len(dresp)) + '.')
+    print('  Number of Dakota response functions = ' + str(len(dresp)))
 
     return dresp
@@ -350 +356 @@
 
     if fline is None or fline == '' or fline.isspace():
-        fline = findline(fidi, 'Moment - based statistics for each response function')
+        fline = findline(fidi, 'Moment-based statistics for each response function')
         return
 
-    print('Reading moment - based statistics for response functions:')
+    print('Reading moment-based statistics for response functions:')
 
     #  skip column headings of moment - based statistics
@@ -376 +382 @@
         dresp[-1].kurtosis = tokens[4]
 
-    print('  Number of Dakota response functions = ' + str(len(dresp)) + '.')
+    print('  Number of Dakota response functions = ' + str(len(dresp)))
 
     return dresp
@@ -433 +439 @@
             dresp[i].stddevci[1, 0] = tokens[4]
 
-    print('  Number of Dakota response functions = ' + str(len(dresp)) + '.')
+    print('  Number of Dakota response functions = ' + str(len(dresp)))
 
     return dresp
@@ -495 +501 @@
                 fline = fidi.readline()
 
-    print('  Number of Dakota response functions = ' + str(len(dresp)) + '.')
+    print('  Number of Dakota response functions = ' + str(len(dresp)))
 
     return dresp
@@ -550 +556 @@
                 fline = fidi.readline()
 
-    print('  Number of Dakota response functions = ' + str(len(dresp)) + '.')
+    print('  Number of Dakota response functions = ' + str(len(dresp)))
 
     return dresp
@@ -565 +571 @@
             return
 
-    print('Reading ' + fline + '.')
+    print('Reading ' + fline)
 
     cmat.title = fline
@@ -657 +663 @@
         #  if importance factors missing, skip to cdf
         if not idvar:
-            print('    Importance Factors not available.')
+            print('    Importance Factors not available')
             dresp[-1].var = []
             dresp[-1].impfac = []
@@ -716 +722 @@
         #  if cdf missing, skip to end of response function
         if not icdf:
-            print('    Cumulative Distribution Function not available.')
+            print('    Cumulative Distribution Function not available')
             dresp[ndresp].cdf = []
             while (fline != '' and not fline.isspace()) and not strncmpi(fline, 'MV Statistics for ', 18) and not strncmp(fline, ' - ', 1):
                 fline = fidi.readline()
 
-    print('  Number of Dakota response functions = ' + str(len(dresp)) + '.')
+    print('  Number of Dakota response functions = ' + str(len(dresp)))
 
     return dresp
@@ -845 +851 @@
         dresp[-1].vum = []
 
-    print('Reading measures for volumetric uniformity.')
+    print('Reading measures for volumetric uniformity')
     fline = fidi.readline()
     fline = fidi.readline()
@@ -881 +887 @@
     [ntokens, tokens] = fltokens(fline)
     method = tokens[2]
-    print('Dakota iterator \'' + str(method) + '\' completed.')
+    print('Dakota iterator \'' + str(method) + '\' completed')
 
     return method
@@ -905 +911 @@
 
     #  issue warning and reset file position
-    print('Warning: findline:str_not_found: String ' + str(string) + ' not found in file.')
+    print('Warning: findline:str_not_found: String ' + str(string) + ' not found in file')
     fidi.seek(ipos, 0)
     return None