main_utils.py

#region "Import of subroutines"
from openpyxl import load_workbook

import os # Miscellaneous operating system interfaces - https://docs.python.org/2.7/library/os.html
import sys # System-specific parameters and functions - https://docs.python.org/2/library/sys.html
from inspect import currentframe, getframeinfo # Get package for monitoring code execution
import pickle # Python object serialization - https://docs.python.org/2/library/pickle.html
from io import StringIO as sio # StringIO — Read and write strings as files - https://docs.python.org/2/library/stringio.html

# Load routines for math and plotting
import matplotlib.pyplot as plt # Matplotlib is a Python 2D plotting library - https://matplotlib.org/
import lmfit
from scipy.interpolate import splprep, splev
from scipy.signal import butter, filtfilt
from scipy.interpolate import interp1d
import numpy as np # NumPy is the fundamental package for scientific computing with Python - https://docs.scipy.org/doc/

import argparse


# ---------------------------------------------------------------------------------------------------------------
# ---- Classes ----
# ---------------------------------------------------------------------------------------------------------------

class ColorData:
    def __init__(self):
        self.Spline = '#000000' #Black
        self.Eng = '#0000FF' #Blue
        self.VolCst = '#008000' #Green
        self.Volkst = '#008000'  # Green
        self.Nu05 = '#8B2E85' #Purple
        self.Tracking = '#FF0000' #red

        self.Orange = '#FF8208'
        self.Turquise = '#23B6B1'
        self.white = '#ffffff'
        self.gray = '#999999'
        self.black = '#000000'
        self.red = '#FF0000'
        self.green = '#00FF00'
        self.blue = '#0000FF'
        self.yellow = '#FFFF00'
        self.cyan = '#00FFFF'
        self.magenta = '#FF00FF'
        self.orange = '#FF9900'
        self.pink = '#FFAAAA'
        self.purple = '#AA00AA'
        self.teal = '#00AA99'
        self.brown = '#996633'
        self.tan = '#FFCC66'
        self.pinkred = '#FF0048'
        self.redorange = '#FF0048'
        self.orangebrown = '#FF8000'
        self.yellowgreen = '#FFFF00'
        self.greencyan = '#00DBA4'
        self.cyanblue = '#00AAFF'
        self.bluemagenta = '#9000FF'
        self.magentapink = '#FF0099'
        #Grayscale
        self.BrightGray = '#EEEEEE'
        self.ChineseSilver = '#CCCCCC'
        self.SpanishGray = '#999999'
        self.GraniteGray = '#666666'
        self.DarkCharcoal = '#333333'

#Generate the subclass for test data sets containing the dataset itself, the data unit and the data description
class TestDataSet:
    def __init__(self, DataSet, DataUnit,DataDescription):
        self.Data = DataSet
        self.Unit = DataUnit
        self.Description = DataDescription

#Class TestPlanData - class containing information from the filled out test-plan spreadsheet for the given test.
class TestPlanData:
    def __init__(self, TestID, GrantaID, TestEngineer, TestDateAndTime, PlanGenerated, PlanLastSaved, PlanLastSavedBy,
                 SpecimenFabricationDate, SpecimenFabricationDateComment, MaterialDesignation,
                 MaterialDesignationComment, NonstandardMaterialDesignation, NonstandardMaterialDesignationComment,
                 MaterialGradeTradeName, MaterialGradeTradeNameComment, ManufacturingProcess,
                 ManufacturingProcessComment, MaterialColor, MaterialColorComment, TotalNumberOfTestSpecimens,
                 TotalNumberOfTestSpecimensComment, TensileTesterProgramme, TensileTesterProgrammeComment, TestType,
                 TestTypeComment, PatternMethod, PatternMethodComment, GeneralComments, Specimendata):

        self.TestID = str(TestID)
        if GrantaID == 'Input Here':
            GrantaID = 0
        self.GrantaID = int(GrantaID)
        self.TestEngineer = str(TestEngineer)
        self.TestDateAndTime = str(TestDateAndTime)
        self.PlanGenerated = str(PlanGenerated)
        self.PlanLastSaved = str(PlanLastSaved)
        self.PlanLastSavedBy = str(PlanLastSavedBy)
        
        self.SpecimenFabricationDate = str(SpecimenFabricationDate)

        if str(SpecimenFabricationDateComment) == 'Input Here':
            self.SpecimenFabricationDateComment = ''
        else:
            self.SpecimenFabricationDateComment = str(SpecimenFabricationDateComment)


        self.MaterialDesignation = str(MaterialDesignation)

        if str(MaterialDesignationComment) == 'Input Here':
            MaterialDesignationComment =''
        else:
            MaterialDesignationComment = str(MaterialDesignationComment)
        self.MaterialDesignationComment = str(MaterialDesignationComment)

        self.NonstandardMaterialDesignation = str(NonstandardMaterialDesignation)

        if str(NonstandardMaterialDesignationComment) == 'Input Here':
            NonstandardMaterialDesignationComment = ''
        else:
            NonstandardMaterialDesignationComment = str(NonstandardMaterialDesignationComment)
        self.NonstandardMaterialDesignationComment = NonstandardMaterialDesignationComment

        self.MaterialGradeTradeName = str(MaterialGradeTradeName)

        if str(MaterialGradeTradeNameComment) == 'Input Here':
            MaterialGradeTradeNameComment = ''
        else:
            MaterialGradeTradeNameComment = str(MaterialGradeTradeNameComment)
        self.MaterialGradeTradeNameComment = MaterialGradeTradeNameComment

        self.ManufacturingProcess = str(ManufacturingProcess)

        if str(ManufacturingProcessComment) == 'Input Here':
            ManufacturingProcessComment = ''
        else:
            ManufacturingProcessComment = str(ManufacturingProcessComment)
        self.ManufacturingProcessComment = ManufacturingProcessComment

        self.MaterialColor = str(MaterialColor)

        if str(MaterialColorComment) == 'Input Here':
            MaterialColorComment = ''
        else:
            MaterialColorComment = str(MaterialColorComment)
        self.MaterialColorComment = MaterialColorComment

        self.TotalNumberOfTestSpecimens = str(TotalNumberOfTestSpecimens)

        if str(TotalNumberOfTestSpecimensComment) == 'Input Here':
            TotalNumberOfTestSpecimensComment = ''
        else:
            TotalNumberOfTestSpecimensComment = str(TotalNumberOfTestSpecimensComment)
        self.TotalNumberOfTestSpecimensComment = TotalNumberOfTestSpecimensComment

        self.TensileTesterProgramme = TensileTesterProgramme

        if str(TensileTesterProgrammeComment) == 'Input Here':
            TensileTesterProgrammeComment = ''
        else:
            TensileTesterProgrammeComment = str(TensileTesterProgrammeComment)
        self.TensileTesterProgrammeComment = TensileTesterProgrammeComment

        if str(TestType) == 'Select':
            TestType = ''
        else:
            TestType = str(TestType)
        self.TestType = TestType

        if str(TestTypeComment) == 'Input Here':
            TestTypeComment = ''
        else:
            TestTypeComment = str(TestTypeComment)
        self.TestTypeComment = TestTypeComment

        if str(GeneralComments) == 'Input Here':
            GeneralComments = ''
        else:
            GeneralComments = str(GeneralComments )

        if str(PatternMethod) == 'Select':
            PatternMethod = ''
        else:
            PatternMethod = str(PatternMethod)
        self.PatternMethod = PatternMethod

        if str(PatternMethodComment) == 'Input Here':
            PatternMethodComment = ''
        else:
            PatternMethodComment = str(PatternMethodComment)
        self.PatternMethodComment = PatternMethodComment

        self.GeneralComments = GeneralComments

        self.Specimendata = {}
        self.Specimendata = Specimendata


# ---------------------------------------------------------------------------------------------------------------
# ---- Subroutines ----
# ---------------------------------------------------------------------------------------------------------------
def getInputVariables():

    ap = argparse.ArgumentParser(description="This code is used to post-process standard DIC test measurements.")
    ap.add_argument(
        'filename', 
        help='The test folder. It is extracted automatically from the command line.')
    ap.add_argument(
        "specimen", 
        type=str, 
        nargs='?',
        default='all',
        help="" \
        "Enter specimen number(s) just after command. " \
        "If multiple samples are entered use space as a separator and enclose in quotes, e.g. '1 2 3'. " \
        "To post-process all specimens, type 'all'. This is also the default value - applied if no value is provided.")

    args = vars(ap.parse_args())
    return args['filename'], args['specimen']


#remove all files with given extension in given folder
def RemoveFiles(Folder, extension):
    test = os.listdir(Folder)
    for item in test:
        if item.endswith(extension):
            os.remove(os.path.join(Folder, item))


#Generate subfolders for the given DIC project.
def GenerateSubfolders(MainFolder):

    if not os.path.isdir(MainFolder + '\\\\Pictures'):
        os.mkdir(MainFolder + '\\\\Pictures')
    else:
        extension = '.png'
        RemoveFiles(MainFolder + '\\\\Pictures', extension)

    if not os.path.isdir(MainFolder + '\\\\PicturePickles'):
        os.mkdir(MainFolder + '\\\\PicturePickles')
    else:
        extension = '.pickle'
        RemoveFiles(MainFolder + '\\\\PicturePickles', extension)

    if not os.path.isdir(MainFolder + '\\\\PostOut'):
        os.mkdir(MainFolder + '\\\\PostOut')
    else:
        extension = '.log'
        RemoveFiles(MainFolder + '\\\\PostOut', extension)
        extension = '.mp4'
        RemoveFiles(MainFolder + '\\\\PostOut', extension)
        extension = '.pickle'
        RemoveFiles(MainFolder + '\\\\PostOut', extension)


#generic routine for fitting a polynomial to data.
def polyfit(x, y, degree):
    results = {}

    coeffs = np.polyfit(x, y, degree)
    results['polynomial'] = coeffs.tolist()

    # r-squared
    p = np.poly1d(coeffs)
    # fit values, and mean
    yhat = p(x)  # or [p(z) for z in x]
    ybar = np.sum(y) / len(y)  # or sum(y)/len(y)
    ssreg = np.sum((yhat - ybar) ** 2)  # or sum([ (yihat - ybar)**2 for yihat in yhat])
    sstot = np.sum((y - ybar) ** 2)  # or sum([ (yi - ybar)**2 for yi in y])
    results['determination'] = ssreg / sstot

    #np.sum((yi - ybar) ** 2)
    return results


def FitSecantModulusFunction(params, x, data):
    '''The FitSecantModulusFunction returns the residuals from the least squares fit +
    the square of any positive slope in the secant fitting function. The secant should be monotonically decreasing
    This is currently obtained by means of introducing a penalty term which penalizes any positive slope
    Currently the parametrization of the secant is: A*strain^2+B*strain+C, with x representing strain'''
    c = params['c']
    b = params['b']
    a = params['a']
    # eps_data = params['eps_data']

    model = a * x * x + b * x + c #secant
    slope = 2 * a * x + b # derivative of secant
    slope[slope < 0.0] = 0.0 # force neg. slopes to zero
    maxSlope = max(slope) # get max of slope
    #FinalResidual = ((data - model) ** 2 * x)
    FinalResidual = ((data - model) ** 2 * x)  #calculate residual
    MaxResidual = max(FinalResidual) #get max of the residuals

    if maxSlope>0.0:
        PenaltyPar = maxSlope / MaxResidual * 1000 #calculate penalty
    else:
        PenaltyPar = 0.0

    FinalResidual[-1] = FinalResidual[-1] * 500.00 #penalize any final values hard

    return FinalResidual + slope * PenaltyPar


def GenerateSecantPolynomial(x, y, polyobject, designation):
    '''Routine for calculating the parameters of the quadratic fitting function of the Secant Modulus
       The routine cals the general minimizer from the lmfit package along with the FitSecantModulusFunction subroutine
       and returns a polyobject with the parameters for the secant modulus.
    '''
    params = lmfit.Parameters() #generate parameters structure

    #add additional parameters: a, b and c
    params.add('c', value=1000, min=0.0)
    params.add('b', value=0.0)  # max=0.0
    params.add('a', value=0.0)

    #minimize residuals and assign result
    minner = lmfit.Minimizer(FitSecantModulusFunction, params, fcn_args=(x, y))
    result = minner.minimize()

    #assign results to local parameters and put these in an array, return the array
    a = result.params['a'].value
    b = result.params['b'].value
    c = result.params['c'].value

    print('GenerateSecantPolynomial: ' + designation + ' done, final chisqr is:' + str(round(result.chisqr, 3)))
    polyobject = np.array([a, b, c])
    return polyobject


#calculate least squares of a one parameter linear regression and return these
def LinearRegression(params, x, data):
    '''The FitSecantModulusFunction returns the residuals from the least squares fit +
    the square of any positive slope in the secant fitting function. The secant should be monotonically decreasing'''
    # a = params['a']
    b = params['b']
    model = b * x
    return ((data - model) ** 2 * x)


#calculate least squares of a two parameter linear regression
def TwoParLinearRegression(params, x, data):
    '''The TwoParLinearRegression returns the residuals from the
    linefit on a set of data'''
    a = params['a']
    b = params['b']
    model = b * x + a
    return ((data - model) ** 2 * x)


def FitPoissonsRatio(Strains, Widthstrains):
    '''calculates estimates of the poisson' ratio based on specimen contraction ratios.'''
    params = lmfit.Parameters()
    params.add('b', value=1.0)  # max=0.0

    # identify optimal poissons ratio in 1.5-3.5% strain interval
    #eliminate any possible doublettes and shrink vectors accordingly
    Shrunk_indices = np.unique(Strains, return_index=True)
    Strains = np.unique(Strains)
    Widthstrains = Widthstrains[Shrunk_indices[1]]

    # calculate poisson's ratios
    entries = len(Widthstrains) - 1
    if entries > 0:
        nuVect = np.zeros(entries + 1)
        Chisqr = np.ones(entries + 1)
        for counter in range(entries + 1 - 2):
            minner = lmfit.Minimizer(LinearRegression, params, fcn_args=(
            abs(Strains[0:entries - counter]), abs(Widthstrains[0:entries - counter])))
            result = minner.minimize()
            if Strains[entries - counter] > 0.017:
                nuVect[counter] = result.params['b'].value
                Chisqr[counter] = result.chisqr
            else:
                nuVect[counter] = 0.4
            
            print("poisson's ratio estimate:" + str(nuVect[counter]) + " result.chisqr" + str(Chisqr[counter]) + " MaxStrain=" + str(
                    Strains[entries - counter]))

        location = np.argmin(Chisqr)
        nu = nuVect[location]
        print("Final Poisson\'s, estimate: " + str(round(nuVect[location], 2)))

    else:
        nu = 0.4
        print("not possible to establish poisson ratio estimate, assuming: 0.4")

    return nu, location


#calculate strain correction due to tensile tester preload.
def FitStrainCorrection(Strains, Force):
    '''Routine for calculating the strain correction due to preload
    this is done using the relation between the force
    and strain slopes.
    '''
    params = lmfit.Parameters()

    params.add('a', value=0.0)
    params.add('b', value=1.0)  # max=0.0

    # strain and force vectors are shrunk in case of non-unique entries.
    Shrunk_indices = np.unique(Strains, return_index=True)
    Strains = np.unique(Strains)
    Force = Force[Shrunk_indices[1]]

    entries = len(Force) - 1
    if entries > 0:
        ABVect = np.zeros(2)
        Chisqr = np.ones(entries + 1)
        #for counter in range(entries + 1 - 2):
        minner = lmfit.Minimizer(TwoParLinearRegression, params, fcn_args=(Force,Strains))
        result = minner.minimize()

        ABVect = [result.params['a'].value, result.params['b'].value]

        StrainCorrection = - ABVect[0]
        ForceStrainSlope = ABVect[1]

    return StrainCorrection, ForceStrainSlope


#Butterworth low pass filter.
def butter_lowpass(cutoff, fs, order=3):
    nyq = 0.5 * fs
    normal_cutoff = cutoff / nyq
    b, a = butter(order, normal_cutoff, btype='low', analog=False)
    return b, a


#Butterworth low pass filter - filtering of data
def butter_lowpass_filtfilt(data, cutoff, fs, order=3):
    b, a = butter_lowpass(cutoff, fs, order=order)
    y = filtfilt(b, a, data)
    return y


# peak-valley detection algorithm
def detect_peaks(x, mph=None, mpd=1, threshold=0, edge='rising',
                 kpsh=False, valley=False, show=False, ax=None):
   
    x = np.atleast_1d(x).astype('float64')
    if x.size < 3:
        return np.array([], dtype=int)
    if valley:
        x = -x
    # find indices of all peaks
    dx = x[1:] - x[:-1]
    # handle NaN's
    indnan = np.where(np.isnan(x))[0]
    if indnan.size:
        x[indnan] = np.inf
        dx[np.where(np.isnan(dx))[0]] = np.inf
    ine, ire, ife = np.array([[], [], []], dtype=int)
    if not edge:
        ine = np.where((np.hstack((dx, 0)) < 0) & (np.hstack((0, dx)) > 0))[0]
    else:
        if edge.lower() in ['rising', 'both']:
            ire = np.where((np.hstack((dx, 0)) <= 0) & (np.hstack((0, dx)) > 0))[0]
        if edge.lower() in ['falling', 'both']:
            ife = np.where((np.hstack((dx, 0)) < 0) & (np.hstack((0, dx)) >= 0))[0]
    ind = np.unique(np.hstack((ine, ire, ife)))
    # handle NaN's
    if ind.size and indnan.size:
        # NaN's and values close to NaN's cannot be peaks
        ind = ind[np.in1d(ind, np.unique(np.hstack((indnan, indnan - 1, indnan + 1))), invert=True)]
    # first and last values of x cannot be peaks
    if ind.size and ind[0] == 0:
        ind = ind[1:]
    if ind.size and ind[-1] == x.size - 1:
        ind = ind[:-1]
    # remove peaks < minimum peak height
    if ind.size and mph is not None:
        ind = ind[x[ind] >= mph]
    # remove peaks - neighbors < threshold
    if ind.size and threshold > 0:
        dx = np.min(np.vstack([x[ind] - x[ind - 1], x[ind] - x[ind + 1]]), axis=0)
        ind = np.delete(ind, np.where(dx < threshold)[0])
    # detect small peaks closer than minimum peak distance
    if ind.size and mpd > 1:
        ind = ind[np.argsort(x[ind])][::-1]  # sort ind by peak height
        idel = np.zeros(ind.size, dtype=bool)
        for i in range(ind.size):
            if not idel[i]:
                # keep peaks with the same height if kpsh is True
                idel = idel | (ind >= ind[i] - mpd) & (ind <= ind[i] + mpd) \
                       & (x[ind[i]] > x[ind] if kpsh else True)
                idel[i] = 0  # Keep current peak
        # remove the small peaks and sort back the indices by their occurrence
        ind = np.sort(ind[~idel])

    if show:
        if indnan.size:
            x[indnan] = np.nan
        if valley:
            x = -x
        #_plot(x, mph, mpd, threshold, edge, valley, ax, ind)

    return ind


#algorithm for finding index closest to the reference value.
def find_idx_nearest_val(array, value):
    # type: (object, object) -> object
    '''Finds the index of value closest to the reference value'''

    idx_sorted = np.argsort(array)
    sorted_array = np.array(array[idx_sorted])
    idx = np.searchsorted(sorted_array, value, side="left")
    if idx >= len(array):
        idx_nearest = idx_sorted[len(array) - 1]
    elif idx == 0:
        idx_nearest = idx_sorted[0]
    else:
        if abs(value - sorted_array[idx - 1]) < abs(value - sorted_array[idx]):
            idx_nearest = idx_sorted[idx - 1]
        else:
            idx_nearest = idx_sorted[idx]
    return idx_nearest


#routing for saving picture in accordance with different modes (pickle format or .png format).
def SaveFigAsPickleAndPngOrInteractive(name, figobject, PicturesDumped, mode):
    '''Function which generat figures as pickle objects + as pictures.
       pictures can be loaded using the below command line
       with open('D:\\ABS_R27_Post\\Specimen4_100mm_min\\Secant_modulus_as_function_of_strain' + '.pickle', 'rb') as handle:
       Picture= pickle.load(handle)'''

    if mode == 'save':
        with open(name[0] + '\\PicturePickles\\' + name[1] + '.pickle', 'wb') as handle:
            pickle.dump(figobject, handle)
        handle.close()

        plt.savefig(name[0] + '\\Pictures\\' + name[1] + '.png')
        PicturesDumped.append(name)
    if mode == 'png':
        plt.savefig(name[0] + '\\Pictures\\' + name[1] + '.png')
        PicturesDumped.append(name)

    return PicturesDumped

#routine for reading Point data exported from VIC-3D, the routine returns a data structure named "Timeseries" which is returned from the routine.
def LoadSpecimenPointData(Mainfolder, subfolderTxtFile, Timeseries, name):

    f = open(os.path.join(Mainfolder, subfolderTxtFile))
    Inputfile = f.read()
    f.close()

    # read the header information into list we can parse - split into words and lower everything
    HeaderInfo = [elem.lower() for elem in Inputfile.splitlines(True)[0:2]]
    NumCategories = HeaderInfo[0].count('"') / 2
    SplittedHeader = HeaderInfo[0].split('"')
    test_delimeter = ","

    if ";" in HeaderInfo[1]:
        test_delimeter = ";"

    # Generate dict with the number of subentries and their names + read these from the SplittedHeader
    NumSubentries = {}
    SubentryNames = {}
    for counter2 in range(0, round(len(SplittedHeader) / 2) - 1):
        NumSubentries[counter2] = SplittedHeader[(counter2 * 2) + 2].count(test_delimeter)
        SubentryNames[counter2] = SplittedHeader[(counter2 * 2) + 1]

    DataCategories = HeaderInfo[1].split(test_delimeter)

    for counter in range(len(DataCategories)):
        DataCategories[counter] = DataCategories[counter].replace('\n', '')

    NumSubentries[0] = NumSubentries[0] - 1
    # get the vector data by using genfromtxt from numpy
    data = np.genfromtxt(sio(Inputfile), delimiter=test_delimeter, skip_header=2)
    # get the transpose of the data
    data = data.transpose()
    # Associate index data

    # Initialize Datacounter
    Datacounter = 1
    # Correct the last number of points in NumSubentries
    NumSubentries[len(NumSubentries) - 1] = NumSubentries[len(NumSubentries) - 1] + 1
    # Loop all the numeric data read and transfer the data into the Timeseries Tree structure.
    for counter2 in range(0, len(NumSubentries)):
        Timeseries[name][SubentryNames[counter2]] = {}
        Timeseries[name][SubentryNames[counter2]]['"index [1]"'] = data[0]
        for counter3 in range(0, NumSubentries[counter2]):
            Timeseries[name][SubentryNames[counter2]][DataCategories[Datacounter]] = data[Datacounter]
            Datacounter = Datacounter + 1

    return Timeseries


#Read the test plan for the given test series
def ReadTestPlan(PlanFolder, Mainfolder, Testplan):

    os.chdir(Mainfolder)
    FileName = Testplan.replace(Testplan[-5:],'').replace('Test_Plan_','')
    wb = load_workbook(os.path.join(PlanFolder, Testplan), read_only=True, data_only=True)
    ws = wb['Test_Plan']

    EndIndex = 17

    PatternMethod = "Select"
    PatternMethodComment = "Input Here"

    # Default values
    TestID = "N/A"
    GrantaID = 0
    TestEngineer = "N/A"
    TestDateAndTime = "N/A"
    PlanGenerated = "N/A"
    PlanLastSaved = "N/A"
    PlanLastSavedBy = "N/A"
    SpecimenFabricationDate = "N/A"
    SpecimenFabricationDateComment = "N/A"
    MaterialDesignation = "N/A"
    MaterialDesignationComment = "N/A"
    NonstandardMaterialDesignation = "N/A"
    NonstandardMaterialDesignationComment = "N/A"
    MaterialGradeTradeName = "N/A"
    MaterialGradeTradeNameComment = "N/A"
    ManufacturingProcess = "N/A"
    ManufacturingProcessComment = "N/A"
    MaterialColor = "N/A"
    MaterialColorComment = "N/A"
    TotalNumberOfTestSpecimens = "N/A"
    TotalNumberOfTestSpecimensComment = "N/A"
    TensileTesterProgramme = "N/A"
    TensileTesterProgrammeComment = "N/A"
    TestType = "N/A"
    TestTypeComment = "N/A"
    PatternMethod = "N/A"
    PatternMethodComment = "N/A"
    GeneralComments = "N/A"


    # Extracting actual data from the test plan
    for counter in range(2, EndIndex):

        if "Test ID" in str(ws['b' + str(counter)].value):
            TestID = ws['d'+str(counter)].value

        if "Granta ID" in str(ws['b' + str(counter)].value):
            GrantaID = ws['d'+str(counter)].value

        if "Test Engineer" in str(ws['b' + str(counter)].value):
            TestEngineer = ws['d'+str(counter)].value

        if "Test date" in str(ws['b' + str(counter)].value):
            TestDateAndTime = ws['d'+str(counter)].value

        if "fabrication date" in str(ws['b' + str(counter)].value):
            SpecimenFabricationDate = ws['d'+str(counter)].value
            SpecimenFabricationDateComment = ws['g'+str(counter)].value

        if "Material Designation" in str(ws['b' + str(counter)].value):
            MaterialDesignation = ws['d'+str(counter)].value
            MaterialDesignationComment = ws['g'+str(counter)].value

        if "Material Grade" in str(ws['b' + str(counter)].value):
            MaterialGradeTradeName = ws['d'+str(counter)].value
            MaterialGradeTradeNameComment = ws['g'+str(counter)].value

        if "Manufacturing process" in str(ws['b' + str(counter)].value):
            ManufacturingProcess = ws['d'+str(counter)].value
            ManufacturingProcessComment = ws['g'+str(counter)].value

        if "Material Color" in str(ws['b' + str(counter)].value):
            MaterialColor = ws['d'+str(counter)].value
            MaterialColorComment = ws['g'+str(counter)].value

        if "Total Number" in str(ws['b' + str(counter)].value):
            TotalNumberOfTestSpecimens = ws['d'+str(counter)].value
            TotalNumberOfTestSpecimensComment = ws['g'+str(counter)].value

        if "Tester Programme" in str(ws['b' + str(counter)].value):
            TensileTesterProgramme = ws['d'+str(counter)].value
            TensileTesterProgrammeComment = ws['g'+str(counter)].value

        if "Test Type" in str(ws['b' + str(counter)].value):
            TestType = ws['d'+str(counter)].value
            TestTypeComment = ws['g'+str(counter)].value

        if "Pattern Method" in str(ws['b' + str(counter)].value):
            PatternMethod = ws['d'+str(counter)].value
            PatternMethodComment = ws['g'+str(counter)].value

    GeneralComments = ws['c'+str(EndIndex+1)].value

    # Extracting the specimen data
    Specimendata = {}

    StartIndex = EndIndex + 3
    EndIndex = StartIndex + int(TotalNumberOfTestSpecimens)

    for RowCount in range(StartIndex, EndIndex, 1):
        Specimendata[RowCount - StartIndex] = {}
        Specimendata[RowCount - StartIndex]["Number"] = str(ws['b' + str(RowCount)].value)
        Specimendata[RowCount - StartIndex]["Specimen Type"] = str(ws['c' + str(RowCount)].value)
        Specimendata[RowCount - StartIndex]["TestSpeed"] = str(ws['d' + str(RowCount)].value)
        Specimendata[RowCount - StartIndex]["AppOrDev."] = str(ws['e' + str(RowCount)].value)
        Specimendata[RowCount - StartIndex]["Orientation"] = str(ws['f' + str(RowCount)].value)
        Specimendata[RowCount - StartIndex]["SpecSpeComments"] = str(ws['g' + str(RowCount)].value)
        Specimendata[RowCount - StartIndex]["Experimental Observations"] = str(ws['h' + str(RowCount)].value)

    Testplaninfo = TestPlanData(
        TestID, GrantaID, TestEngineer, TestDateAndTime, PlanGenerated, PlanLastSaved, PlanLastSavedBy,
        SpecimenFabricationDate, SpecimenFabricationDateComment, MaterialDesignation,
        MaterialDesignationComment, NonstandardMaterialDesignation,
        NonstandardMaterialDesignationComment, MaterialGradeTradeName,
        MaterialGradeTradeNameComment, ManufacturingProcess, ManufacturingProcessComment,
        MaterialColor, MaterialColorComment, TotalNumberOfTestSpecimens,
        TotalNumberOfTestSpecimensComment, TensileTesterProgramme,
        TensileTesterProgrammeComment, TestType, TestTypeComment, PatternMethod, PatternMethodComment,
        GeneralComments, Specimendata)

    return Testplaninfo, FileName


#This routine reads the Ascii file from the TO tensile testing machine and returns the data in reals + vectors
def LoadspecimenTOdata(Origin, Mainfolder, TOfile, InitialThickness, InitialWidth):

    f = open(os.path.join(Origin, Mainfolder, TOfile))
    Inputfile = f.read().replace(',', '')
    f.close()

    HeaderInfo = [elem.lower() for elem in Inputfile.splitlines(True)[0:2]]
    NumCategories = HeaderInfo[0].count(';')
    SplittedHeader = HeaderInfo[0].split(';')
    SplittedHeaderData = HeaderInfo[1].split(';')
   
    for counter in range(0, NumCategories):
        if 'width' in SplittedHeader[counter].lower():
            InitialWidth = float(SplittedHeaderData[counter])
        if 'thickness' in SplittedHeader[counter].lower():
            InitialThickness = float(SplittedHeaderData[counter])

    HeaderInfo2 = [elem.lower() for elem in Inputfile.splitlines(True)[2:3]][0].split(';')
    inputarray = np.genfromtxt(sio(Inputfile), delimiter=';', skip_header=3, max_rows=len(Inputfile.splitlines(True)))
    Transpose = np.transpose(inputarray)

    HorizonTime = Transpose[[i for i, elem in enumerate(HeaderInfo2) if 'time' in elem][0]]-Transpose[[i for i, elem in enumerate(HeaderInfo2) if 'time' in elem][0]][0]
    HorizonForce = Transpose[[i for i, elem in enumerate(HeaderInfo2) if 'force' in elem][0]]
    HorizonDisplacement = Transpose[[i for i, elem in enumerate(HeaderInfo2) if 'position' in elem][0]]

    if len(HorizonTime) > 1048000:
        HorizonTime = HorizonTime[0::5]
        HorizonForce = HorizonForce[0::5]
        HorizonDisplacement = HorizonDisplacement[0::5]

    return InitialWidth, InitialThickness, HorizonTime, HorizonForce, HorizonDisplacement


#algorithm for replacing NAN values with interpolated data instead.
def fill_nans_scipy1(padata, pkind='linear'):
    """
    Interpolates data to fill nan values
    Parameters:
        padata : nd array
            source data with np.NaN values

    Returns:
        nd array
            resulting data with interpolated values instead of nans

    Requires: from scipy.interpolate import interp1d
    """
    aindexes = np.arange(padata.shape[0])
    agood_indexes, = np.where(np.isfinite(padata))
    f = interp1d(agood_indexes
                 , padata[agood_indexes]
                 , bounds_error=False
                 , copy=False
                 , fill_value="extrapolate"
                 , kind=pkind)
    return f(aindexes)


#routine which identifies the ID's of the end of the report
def Identify_Image_indixes(ForcePeaks, startofRange, endofActiveRange, Lineseries, Timeseries, Mainfolder):

    filelist_old = os.listdir(Mainfolder)
    filelist_old = [x.lower() for x in filelist_old]

    filelist = []
    for flist in filelist_old: # filelist[:] makes a copy of filelist.
        if flist.endswith(".tif"):
            new_flist = flist.split("-")
            if new_flist[1] != "flir0" and new_flist[1] != "fcal0":
                filelist.append(flist)

    #limit = endofActiveRange*2 + 2
    #filelist = filelist[startofRange:limit]

    Imagedigits=len(Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][0][0])

    ImageInformation = {}
    ImageInformation[0] = {}
    ImageInformation[0]['index'] = 0
    ImageInformation[0]['OrgIndex'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][0][0]
    ImageInformation[0]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][0][1]
    ImageInformation[0]['Force'] = str(round(float(Timeseries[list(Timeseries.keys())[0]]['Post']['Gen']['Force'].Data[0]),2))
    ImageInformation[0]['dispTester'] = str(round(float(Timeseries[list(Timeseries.keys())[0]]['Post']['Gen']['dispTester'].Data[0]),3))
    ImageInformation[0]['Neckstrain'] = str(float(round(Timeseries[list(Timeseries.keys())[0]]['Post']['Volkst']['Strain_VolConstant'].Data[0],2)))

    if Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + ImageInformation[0]['OrgIndex'] + '_' +\
            ImageInformation[0]['OrgCamera'] + '.tif' in filelist:
        ImageInformation[0]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + \
                                       ImageInformation[0]['OrgIndex'] + '_' + ImageInformation[0]['OrgCamera'] + '.tif'
    else:
        print("DIC Image for report not found", "a picture searched for in he DIC folder was not found!! - ending program, the path is: "+
                               Mainfolder+Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + ImageInformation[0]['OrgIndex'] + '_' +
                               ImageInformation[0]['OrgCamera'] + '.tif')

    for counter in range(1, len(ForcePeaks)+1):
        ImageInformation[counter] = {}
        ImageInformation[counter]['index'] = ForcePeaks[counter-1]
        ImageInformation[counter]['OrgIndex'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][ForcePeaks[counter-1]][0]
        ImageInformation[counter]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][ForcePeaks[counter-1]][1]
        ImageInformation[counter]['Force'] = str(round(float(Timeseries[list(Timeseries.keys())[0]]['Post']['Gen']['Force'].Data[ForcePeaks[counter-1]]),2))
        ImageInformation[counter]['dispTester'] = str(round(float(Timeseries[list(Timeseries.keys())[0]]['Post']['Gen']['dispTester'].Data[ForcePeaks[counter-1]]),3))
        ImageInformation[counter]['Neckstrain'] = str(float(round(Timeseries[list(Timeseries.keys())[0]]['Post']['Volkst']['Strain_VolConstant'].Data[ForcePeaks[counter-1]], 2)))

        if Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + ImageInformation[counter]['OrgIndex'] + '_' +ImageInformation[counter]['OrgCamera'] + '.tif' in filelist:
            ImageInformation[counter]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + ImageInformation[counter]['OrgIndex'] + '_' + \
                                                 ImageInformation[counter]['OrgCamera'] + '.tif'
        else:
            print("DIC Image for report not found","a picture searched for in he DIC folder was not found!! "
                                                                    "- ending program, the path is: " + Mainfolder)

    if endofActiveRange != ImageInformation[counter]['index']:
        ImageInformation[counter+1] = {}
        ImageInformation[counter+1]['index'] = endofActiveRange
        ImageInformation[counter+1]['OrgIndex'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][0]
        ImageInformation[counter+1]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
        ImageInformation[counter+1]['Force'] = str(round(float(Timeseries[list(Timeseries.keys())[0]]['Post']['Gen']['Force'].Data[endofActiveRange]),2))
        ImageInformation[counter+1]['dispTester'] = str(round(float(Timeseries[list(Timeseries.keys())[0]]['Post']['Gen']['dispTester'].Data[endofActiveRange]),3))
        ImageInformation[counter+1]['Neckstrain'] = str(float(round(Timeseries[list(Timeseries.keys())[0]]['Post']['Volkst']['Strain_VolConstant'].Data[endofActiveRange], 2)))

        if Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + ImageInformation[counter+1]['OrgIndex'] + '_' +ImageInformation[counter+1]['OrgCamera'] + '.tif' in filelist:
            ImageInformation[counter+1]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' +ImageInformation[counter+1]['OrgIndex'].zfill(Imagedigits) + '_' + \
                                                     ImageInformation[counter+1]['OrgCamera'] + '.tif'
        else:
            print("DIC Image for report not found","a picture searched for in he DIC folder was not found!! - ending program, the path is: " + Mainfolder)

        counter = counter + 2
    else:
        counter = counter + 1

    if ImageInformation[0]['OrgCamera']=='0':
        Check = '1'
    else:
        Check = '0'

    removelist =[]
    for counter2 in range(0,len(filelist),1):
        if '_'+Check+'.tif' in filelist[counter2]:
            removelist.append(filelist[counter2])

    filelist = list(set(filelist) - set(removelist))
    indices = np.zeros(len(filelist)).astype(int)

    for counter2 in range(0, len(filelist), 1):
        #print(filelist[counter2])
        indices[counter2]=int(filelist[counter2].replace(Lineseries[list(Lineseries.keys())[0]]['VicProjectname']+'-', '').replace('_'+ImageInformation[0]['OrgCamera']+'.tif', ''))
    MissingPics = sorted(i for i in indices if i > endofActiveRange)

    if len(MissingPics)>0:
        if len(MissingPics)==1:
            ImageInformation[counter] = {}
            ImageInformation[counter]['index'] = 'N.A.'
            ImageInformation[counter]['OrgIndex'] = str(MissingPics[0])
            ImageInformation[counter]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
            ImageInformation[counter]['Force'] = 'N.A.'
            ImageInformation[counter]['dispTester'] = 'N.A.'
            ImageInformation[counter]['Neckstrain'] = 'N.A.'
            ImageInformation[counter]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' +str(MissingPics[0]).zfill(Imagedigits) + '_' + ImageInformation[0]['OrgCamera'] + '.tif'
        else:
            if len(MissingPics)==2:
                ImageInformation[counter] = {}
                ImageInformation[counter]['index'] = 'N.A.'
                ImageInformation[counter]['OrgIndex'] = str(MissingPics[0])
                ImageInformation[counter]['OrgCamera'] = \
                Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
                ImageInformation[counter]['Force'] = 'N.A.'
                ImageInformation[counter]['dispTester'] = 'N.A.'
                ImageInformation[counter]['Neckstrain'] = 'N.A.'
                ImageInformation[counter]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + str(MissingPics[0]).zfill(Imagedigits) + '_' + ImageInformation[counter]['OrgCamera'] + '.tif'
                
                ImageInformation[counter+1] = {}
                ImageInformation[counter+1]['index'] = 'N.A.'
                ImageInformation[counter+1]['OrgIndex'] = str(MissingPics[1])
                ImageInformation[counter+1]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
                ImageInformation[counter+1]['Force'] = 'N.A.'
                ImageInformation[counter+1]['dispTester'] = 'N.A.'
                ImageInformation[counter+1]['Neckstrain'] = 'N.A.'
                ImageInformation[counter+1]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + str(MissingPics[1]).zfill(Imagedigits) + '_' + ImageInformation[counter+1]['OrgCamera'] + '.tif'
            else:
                ImageInformation[counter] = {}
                ImageInformation[counter]['index'] = 'N.A.'
                ImageInformation[counter]['OrgIndex'] = str(MissingPics[0])
                ImageInformation[counter]['OrgCamera'] = \
                Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
                ImageInformation[counter]['Force'] = 'N.A.'
                ImageInformation[counter]['dispTester'] = 'N.A.'
                ImageInformation[counter]['Neckstrain'] = 'N.A.'
                ImageInformation[counter]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + str(MissingPics[0]).zfill(Imagedigits) + '_' + ImageInformation[counter]['OrgCamera'] + '.tif'
                
                ImageInformation[counter+1] = {}
                ImageInformation[counter+1]['index'] = 'N.A.'
                ImageInformation[counter+1]['OrgIndex'] = str(MissingPics[1])
                ImageInformation[counter+1]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
                ImageInformation[counter+1]['Force'] = 'N.A.'
                ImageInformation[counter+1]['dispTester'] = 'N.A.'
                ImageInformation[counter+1]['Neckstrain'] = 'N.A.'
                ImageInformation[counter+1]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + str(MissingPics[1]).zfill(Imagedigits) + '_' + ImageInformation[counter+1]['OrgCamera'] + '.tif'
                
                ImageInformation[counter+2] = {}
                ImageInformation[counter+2]['index'] = 'N.A.'
                ImageInformation[counter+2]['OrgIndex'] = str(MissingPics[-1])
                ImageInformation[counter+2]['OrgCamera'] = Lineseries[list(Lineseries.keys())[0]]['OrgImageArray'][endofActiveRange][1]
                ImageInformation[counter+2]['Force'] = 'N.A.'
                ImageInformation[counter+2]['dispTester'] = 'N.A.'
                ImageInformation[counter+2]['Neckstrain'] = 'N.A.'
                ImageInformation[counter+2]['image'] = Mainfolder + Lineseries[list(Lineseries.keys())[0]]['VicProjectname'] + '-' + str(MissingPics[-1]).zfill(Imagedigits) + '_' + ImageInformation[counter+2]['OrgCamera'] + '.tif'

    return ImageInformation