Merge branch 'develop' of git.paveit.de:paveit/lib-paveit into develop

src/paveit/labtest/__init__.py

@@ -1,5 +1,8 @@
 from .base import DataSineLoad
+from .citt import *
 from .citt import CITTBase
+from .citt_fatigue import *
+from .dsv import *
 
 __all__ = ['DataSineLoad',
             'CITTBase'

src/paveit/labtest/base.py

@@ -4,6 +4,7 @@ import logging
 
 import numpy as np
 import pandas as pd
+
 from paveit.analysis import fit_cos
 from paveit.functions import calc_nu
 from paveit.helper import calc_hash_of_bytes, get_minio_client_processing
@@ -224,6 +225,7 @@ class DataSineLoad():
 
         for name in self.val_header_names:
             if not name in self.metadata:
+                self._logger.error(f'{name} not found')
                 raise
 
     def _post_string_to_float(self):

src/paveit/labtest/citt.py

@@ -10,8 +10,96 @@ from paveit import calc_nu, fit_cos
 from paveit.datamodels import CITTSiffness, CITTSiffnessResults
 from paveit.io import read_geosys
 from paveit.labtest import DataSineLoad
+from paveit.labtest.citt_fatigue import CittAnalyseFatigue
 
 
+def calc_E(data, metadata, columns_analyse):
+
+    data.index = data.index - data.index[0]
+    
+    res_temp = {}
+
+    x = data.index.values
+
+    freq = np.round(float(data['f'].unique()), 2)
+    sigma = float(data['sigma'].unique())
+    temperature = float(data['T'].unique())
+    
+    for idxcol, col in enumerate(columns_analyse):
+
+        if not col in data.columns: continue
+
+        y = data[col].values
+        res = fit_cos(x, y, freq=freq)
+
+        for key, value in res.items():
+            res_temp[f'fit_{col}_{key}'] = value
+            
+
+        # analyse cycle data
+        
+        cycle_min = []
+        cycle_max = []
+        cycle_mean = []
+        cycle_diff = []
+        
+        for N, data_cycle in data.groupby('N'):
+            y = data_cycle[col].values
+            
+            cycle_min.append(y.min())
+            cycle_max.append(y.max())
+            cycle_mean.append(y.mean())
+            cycle_diff.append(cycle_max[-1] - cycle_min[-1])
+            
+        res_temp[f'fit_{col}_cycle_min'] = cycle_min
+        res_temp[f'fit_{col}_min'] = np.mean(cycle_min)
+        res_temp[f'fit_{col}_min_std'] = np.std(cycle_min)
+        res_temp[f'fit_{col}_min_diff_rel'] = (np.max(cycle_min) - np.min(cycle_min))/np.mean(cycle_min)
+        res_temp[f'fit_{col}_cycle_max'] = cycle_max
+        res_temp[f'fit_{col}_max'] = np.mean(cycle_max)
+        res_temp[f'fit_{col}_max_std'] = np.std(cycle_max)
+        res_temp[f'fit_{col}_max_diff_rel'] = (np.max(cycle_max) - np.min(cycle_max))/np.mean(cycle_max)
+        res_temp[f'fit_{col}_cycle_mean'] = cycle_mean
+        res_temp[f'fit_{col}_mean'] = np.mean(cycle_mean)
+        res_temp[f'fit_{col}_mean_std'] = np.std(cycle_mean)
+        res_temp[f'fit_{col}_mean_diff_rel'] = (np.max(cycle_mean) - np.min(cycle_mean))/np.mean(cycle_mean)
+        res_temp[f'fit_{col}_cycle_diff'] = cycle_diff
+        res_temp[f'fit_{col}_diff'] = np.mean(cycle_diff)
+        res_temp[f'fit_{col}_diff_std'] = np.std(cycle_diff)
+        res_temp[f'fit_{col}_diff_diff_rel'] = (np.max(cycle_diff) - np.min(cycle_diff))/np.mean(cycle_diff)
+
+    # add more metadata
+    res_temp['f_set'] = freq
+    res_temp['sigma_set'] = sigma
+    res_temp['T_set'] = temperature
+
+    res_temp['N_from'] = data['N'].min()
+    res_temp['N_to'] = data['N'].max()
+
+    res_temp['n_samples_per_cycle'] = int(
+        len(data) / (res_temp['N_to'] - res_temp['N_from'] + 1))
+
+    ## Stiffness
+    deltaF = res_temp['fit_F_amp']
+    deltaU = res_temp['fit_s_hor_sum_amp']
+    
+    h = float(metadata['speciment_height'])
+    d = float(metadata['speciment_diameter'])
+    
+    nu = calc_nu(temperature)
+    res_temp['nu'] = nu
+    
+    #nach TP Asphalt 26
+    res_temp['stiffness'] = deltaF /(h * deltaU) * (4.0/np.pi -1 + nu)
+    
+    ## Elastische hori. Dehnung                
+    res_temp['el_strains'] = 2*2*deltaU/d * (1+3*nu)/(4 + np.pi*nu - np.pi) * 1000.0 # 2*2 daher, da deltaU nur Ampl. nicht Gesamtkraft ist
+    
+    # TODO: Überarbeiten und erweitern (ISSUE #2)
+    res_temp['phase'] = res_temp['fit_F_phase'] - res_temp['fit_s_hor_sum_phase']
+    
+    return res_temp
+
 class CITTBase(DataSineLoad):
     
     def _set_parameter(self):
@@ -837,6 +925,7 @@ class CITT_LaborHart(CITTBase):
             'speciment_diameter': ['Probendurchmesser'],
             'speciment_height': ['Probenhöhe'],
             'speciment_name': ['Probenbezeichnung'],
+            'f': ['FREQUENZ_x10'],
         }  #list of names
 
         self.data_column_names = {
@@ -846,7 +935,7 @@ class CITT_LaborHart(CITTBase):
             's_hor_1': ['SENSOR_EXT_µm_x10'],
             's_hor_2': ['SENSOR_S4_µm_x10'],
             's_piston': ['POSITION_µm'],
-            'N': ['Zyklen_fortlaufend'],
+            'N': ['Zyklen_fortlaufend', 'Zyklen'],
         }
 
     def _process_data(self):
@@ -918,6 +1007,9 @@ class CITT_LaborHart(CITTBase):
         self._logger.info(self.metadata)
         self._logger.info(self.data.head())
         
+class CITT_fatigue_LaborHart(CittAnalyseFatigue ,CITT_LaborHart):
+    pass
+
 
 class CITT_BAGKoeln(CITTBase):
 

src/paveit/labtest/citt_fatigue.py

@@ -0,0 +1,91 @@
+import numpy as np
+import pandas as pd
+
+from paveit.labtest.citt import calc_E
+
+
+class CittAnalyseFatigue():
+    
+    def _fit_split_data(self):
+        
+        data_exp = []
+        
+        N = self.data['N'].unique()
+        N = np.array(N)
+        gaps = N[1:][np.diff(N)>1]
+
+        for i,gap in enumerate(gaps):
+            print(i, gap)
+            if i == 0:
+                f = self.data['N']<gap
+            elif i == len(gaps):
+                f = self.data['N']>=gap
+            else:
+                f = (self.data['N']>=gaps[i-1]) & (self.data['N']<gap)
+        
+            # filter data by geps
+            d = self.data[f]
+            
+            # get 5 cycles
+            if i == 0:
+                f = (d['N']>=98) & (d['N']<=102)
+            else:
+                Nsel = d['N'].unique()
+                f = (d['N']>=Nsel[-5]) & (d['N']<=Nsel[-1])
+            
+            d = d[f]
+            
+            
+            
+            data_exp.append(d)
+        
+        
+        self.data = data_exp
+        
+    
+    def _fit_select_data(self):
+        
+        ''' analyse data
+
+        '''
+        pass 
+        
+    def _calc(self):
+        print('calc fatigue')
+        
+        print(self.metadata)
+        
+        fit = []
+        
+        # Je Aufzeichnungsintervall
+        
+        
+        for i, d in enumerate(self.data):
+            
+            try:
+                res = calc_E(d, metadata=self.metadata, columns_analyse=['F', 's_hor_sum'])
+            
+                res['idx'] = i
+                
+                res['energy_ratio'] = res['stiffness']*np.round(res['N_from'] + (res['N_to'] - res['N_from'])/2, 0)
+            
+                fit.append(res)
+            except:
+                raise
+            
+        self.fit_single_results = pd.DataFrame.from_records(fit)
+        
+        EN_max = self.fit_single_results['energy_ratio'].max()
+        
+        sel_f = self.fit_single_results[(self.fit_single_results['energy_ratio']>=0.8*EN_max) & (self.fit_single_results['energy_ratio']<=1.2*EN_max)]
+        par = np.polyfit(sel_f['N_from'], sel_f['energy_ratio'], 4)
+        
+        x = np.arange(sel_f['N_from'].min(),sel_f['N_from'].max(), 1)
+        y = np.polyval(par, x)
+
+        Nmakro = x[y.argmax()]
+
+        self.fit = {'Nmakro': Nmakro, 
+                        'energy_ratio_max': y.max(),
+                        'par_fit': par,
+                        'epislon_elast_98': self.fit_single_results.iloc[0]['el_strains']}

src/paveit/labtest/dsv.py

@@ -0,0 +1,394 @@
+import io
+import os
+from csv import reader
+
+import numpy as np
+import pandas as pd
+from bson import ObjectId
+
+from paveit import calc_nu, fit_cos
+from paveit.datamodels import CITTSiffness, CITTSiffnessResults
+from paveit.io import read_geosys
+from paveit.labtest import DataSineLoad
+
+
+class TP25A1base(DataSineLoad):
+    
+    def _set_parameter(self):
+        self._logger.debug('run _set_parameter')
+
+        self.split_data_based_on_parameter = ['T', 'sigma', 'f']
+
+        self.col_as_int = ['N']
+        
+        self.col_as_float = [
+            'T', 'F', 's_piston', 's_hor_1', 'f', 's_hor_1', 's_hor_2'
+        ]
+
+        self.val_col_names = [
+            'time', 'T', 'f', 'sigma', 'N', 'F', 's_hor_1', 's_hor_2', 's_piston'
+        ]
+
+        self.columns_analyse = [
+            'F', 's_hor_sum', 's_hor_1', 's_hor_2', 's_piston'
+        ]
+
+        self.round_values = [('T', 3), ('sigma', 3)]
+
+        # Header names after standardization; check if exists
+        self.val_header_names = ['speciment_height', 'speciment_diameter']
+
+        self.number_of_load_cycles_for_analysis = 5
+
+        #Dummy Data, replace in Machine Config
+        self.meta_names_of_parameter = {
+            'sigma': ['Max. Spannung']
+        }  #list of names
+
+        #Dummy Data, replace in Machine Config
+        self.data_column_names = {
+            'time': ['Time Series'],
+            'F': ['Load Series'],
+            's_hor_1': ['LVDT1 Series'],
+            's_hor_2': ['LVDT2 Series'],
+        }
+
+    def _sel_df(self, df, num=5, shift=-1):
+        
+        print(df.head())
+
+        N = df['N'].unique()
+        n_N = len(N)
+        max_N = max(N)
+        min_N = min(N)
+
+        freq = float(df['f'].unique()[0])
+        
+        # define cycles to select
+        if freq == 10.0:
+            Nfrom = 98
+            Nto = 103
+        elif freq == 5.0:
+            Nfrom = 93
+            Nto = 97
+        elif freq == 3.0:
+            Nfrom = 43
+            Nto = 47
+        elif freq == 1.0:
+            Nfrom = 13
+            Nto = 17
+        elif freq == 0.3:
+            Nfrom = 8
+            Nto = 12
+        elif freq == 0.1:
+            Nfrom = 3
+            Nto = 7
+        else:
+            Nfrom = None
+            Nto = None
+            
+            
+        self._logger.debug(f'{min_N}, {max_N}, {n_N}, {num}, {shift}')
+        self._logger.debug(f'Frequenz: {freq}, Nfrom: {Nfrom}, Nto: {Nto}')
+
+        # Fall 1: nur num Lastwechsel
+        if n_N < num:
+            df_sel = None
+        elif n_N == num:
+            df_sel = df
+        # Fall 2: nicht alle LW in Datei
+        elif (max_N < Nto) & (n_N > num):
+
+            df_sel = df[(df['N'] >= N[-num + shift])
+                        & (df['N'] <= N[-1 + shift])]
+
+        # Fall 3: Auswahl wie oben definiert
+        elif (Nfrom >= min_N) & (Nto < max_N):
+            df_sel = df[(df['N'] >= Nfrom) & (df['N'] <= Nto)]
+        # Fall 4: Auswahl unbekannt
+        else:
+            df_sel = None
+
+        return df_sel
+
+    def _fit_select_data(self):
+        """ 
+        select N load cycles from original data 
+        (a): Based on window of TP-Asphalt
+        (b) last N cycles
+        
+        """
+
+        self._logger.debug('run _fit_select_data')
+
+        self.max_N_in_data = []
+
+        if not isinstance(self.data, list):
+            if self.number_of_load_cycles_for_analysis > 1:
+
+                self.max_N_in_data.append(self.data['N'].max())
+
+                df_sel = [
+                    self._sel_df(self.data,
+                                 num=self.number_of_load_cycles_for_analysis)
+                ]
+            else:
+                df_sel = [self.data]
+
+        else:
+            df_sel = []
+            for d in self.data:
+
+                self.max_N_in_data.append(d['N'].max())
+
+                if self.number_of_load_cycles_for_analysis > 1:
+                    d_sel = self._sel_df(
+                        d, num=self.number_of_load_cycles_for_analysis)
+                else:
+                    d_sel = d
+
+                df_sel.append(d_sel)
+
+        # replace data
+        self.data = df_sel
+
+    def _calc(self):
+
+        self._logger.info('run _calc CITT')
+        print('run CITT')
+
+        self.fit = []
+        for idx_data, data in enumerate(self.data):
+
+            if data is None: continue
+            if len(data) < 10: continue
+            
+            try:
+                
+                self._logger.debug(f'run fit on subset {idx_data}')
+
+                data.index = data.index - data.index[0]
+
+                res_temp = {}
+                res_temp['idx'] = idx_data
+
+                x = data.index.values
+
+                freq = np.round(float(data['f'].unique()), 2)
+                sigma = float(data['sigma'].unique())
+                temperature = float(data['T'].unique())
+
+                for idxcol, col in enumerate(self.columns_analyse):
+
+                    if not col in data.columns: continue
+
+                    y = data[col].values
+                    res = fit_cos(x, y, freq=freq)
+
+                    for key, value in res.items():
+                        res_temp[f'fit_{col}_{key}'] = value
+                        
+
+                    # analyse cycle data
+                    
+                    cycle_min = []
+                    cycle_max = []
+                    cycle_mean = []
+                    cycle_diff = []
+                    
+                    for N, data_cycle in data.groupby('N'):
+                        y = data_cycle[col].values
+                        
+                        cycle_min.append(y.min())
+                        cycle_max.append(y.max())
+                        cycle_mean.append(y.mean())
+                        cycle_diff.append(cycle_max[-1] - cycle_min[-1])
+                        
+                    res_temp[f'fit_{col}_cycle_min'] = cycle_min
+                    res_temp[f'fit_{col}_min'] = np.mean(cycle_min)
+                    res_temp[f'fit_{col}_min_std'] = np.std(cycle_min)
+                    res_temp[f'fit_{col}_min_diff_rel'] = (np.max(cycle_min) - np.min(cycle_min))/np.mean(cycle_min)
+                    res_temp[f'fit_{col}_cycle_max'] = cycle_max
+                    res_temp[f'fit_{col}_max'] = np.mean(cycle_max)
+                    res_temp[f'fit_{col}_max_std'] = np.std(cycle_max)
+                    res_temp[f'fit_{col}_max_diff_rel'] = (np.max(cycle_max) - np.min(cycle_max))/np.mean(cycle_max)
+                    res_temp[f'fit_{col}_cycle_mean'] = cycle_mean
+                    res_temp[f'fit_{col}_mean'] = np.mean(cycle_mean)
+                    res_temp[f'fit_{col}_mean_std'] = np.std(cycle_mean)
+                    res_temp[f'fit_{col}_mean_diff_rel'] = (np.max(cycle_mean) - np.min(cycle_mean))/np.mean(cycle_mean)
+                    res_temp[f'fit_{col}_cycle_diff'] = cycle_diff
+                    res_temp[f'fit_{col}_diff'] = np.mean(cycle_diff)
+                    res_temp[f'fit_{col}_diff_std'] = np.std(cycle_diff)
+                    res_temp[f'fit_{col}_diff_diff_rel'] = (np.max(cycle_diff) - np.min(cycle_diff))/np.mean(cycle_diff)
+
+                # add more metadata
+                res_temp['f_set'] = freq
+                res_temp['sigma_set'] = sigma
+                res_temp['T_set'] = temperature
+
+                res_temp['N_from'] = data['N'].min()
+                res_temp['N_to'] = data['N'].max()
+                res_temp['N_tot'] = self.max_N_in_data[idx_data]
+
+                res_temp['n_samples_per_cycle'] = int(
+                    len(data) / (res_temp['N_to'] - res_temp['N_from'] + 1))
+
+                ## Stiffness
+                deltaF = res_temp['fit_F_amp']
+                deltaU = res_temp['fit_s_hor_sum_amp']
+                
+                h = float(self.metadata['speciment_height'])
+                d = float(self.metadata['speciment_diameter'])
+                
+                nu = calc_nu(temperature)
+                res_temp['nu'] = nu
+                
+                print(deltaF, deltaU, h, d, nu, np.pi)
+
+                #nach TP Asphalt 26
+                res_temp['stiffness'] = deltaF /(h * deltaU) * (4.0/np.pi -1 + nu)
+                
+                ## Elastische hori. Dehnung                
+                res_temp['el_strains'] = 2*2*deltaU/d * (1+3*nu)/(4 + np.pi*nu - np.pi) * 1000.0 # 2*2 daher, da deltaU nur Ampl. nicht Gesamtkraft ist
+                
+                # TODO: Überarbeiten und erweitern (ISSUE #2)
+                res_temp['phase'] = res_temp['fit_F_phase'] - res_temp['fit_s_hor_sum_phase']
+                
+            except Exception as e:
+                self._logger.exception(e)
+                res_temp = None
+                
+            self._logger.debug(res_temp)
+
+            self.fit.append(res_temp)
+
+        self.fit = pd.DataFrame.from_records(self.fit)
+
+        self.fit = self.fit.reset_index(drop=True).set_index('idx')
+
+        #self.fit = self.fit.set_index(['T', 'f', 'sigma'])
+
+        nsamples = len(self.fit)
+        self._logger.info(f'fitting finished, add {nsamples} samples')
+        self._logger.debug(self.fit['stiffness'])
+
+    def save(self,
+             task_id: ObjectId,
+             meta: dict = {}
+             ):
+        """
+        save results to mongodb
+        """
+
+        if not hasattr(self, 'fit'):
+            raise
+
+        # precheck data and results
+        #assert len(self.data) == len(self.fit)
+
+        for idx_fit, fit in self.fit.iterrows():
+            data = self.data[idx_fit]
+
+            meta['filehash'] = self.filehash
+            meta['task_id'] = task_id
+            
+            if not self.metadata['speciment_name'] == None:
+                meta['speciment_name'] = self.metadata['speciment_name']
+            else:
+                meta['speciment_name'] = self.filename
+            
+            meta['speciment_diameter'] = self.metadata['speciment_diameter']
+            meta['speciment_height'] = self.metadata['speciment_height']
+            
+
+            #check if result in db
+            #n = CITTSiffness.objects(**meta).count()
+            #print(n)
+
+            # write data
+            data_dict = fit.to_dict()
+            data_dict.update(meta)
+
+            # remove 'fit_' from keys:
+            for key in list(data_dict.keys()):
+                if key.startswith('fit_'):
+                    data_dict[key[4:]] = data_dict[key]
+                    data_dict.pop(key)
+
+            # rename fields
+
+            def rename_field(d, old, new):
+                d[new] = d[old]
+                d.pop(old)
+
+            f = CITTSiffnessResults(**data_dict).save()
+
+            # required data
+            data_out = dict(
+                time=data.index,
+                F=list(data['F']),
+                N=list(data['N']),
+                s_hor_1=list(data['s_hor_1']),
+                s_hor_2=list(data['s_hor_2']),
+                s_hor_sum=list(data['s_hor_sum']),
+            )
+            
+            self._logger.debug(f'columns data,  {data.columns}')
+
+            # add optional datas
+            for col in ['s_piston']:
+                if col in data.columns:
+                    self._logger.debug(f'add {col} to output data')
+                    data_out[col] = list(data[col])
+
+            outkeys = list(data_out.keys())
+            self._logger.debug(f'write raw data to db,  {outkeys}')
+            g = CITTSiffness(result=f.id, **data_out).save()
+
+class TP25A1_TUDresdenWille(TP25A1base):
+
+    def _define_units(self):
+
+        self.unit_s = 1 / 1000.  #mm
+        self.unit_F = 1.0  #N
+        self.unit_t = 1.  #s
+
+    def update_parameter(self):
+
+        self.meta_names_of_parameter = {
+            'sigma': ['Oberspannung'],
+            'sigma_min': ['Unterspannung'],
+            'f': ['Frequenz'],
+            'T': ['Versuchstemperatur'],
+            't': ['Zeit'],
+            't_pulse': ['Impulsdauer'],
+            't_break': ['Lastpauseit'],
+            'speciment_diameter': ['PK-Durchmesser'],
+            'speciment_height': ['PK-Höhe'],
+            'punch_diameter': ['Stempeldurchmesser'],
+            'speciment_name': ['Probekörperbezeichnung'],
+        }  #list of names
+
+        self.data_column_names = {
+            'time': ['Zeit'],
+            'F': ['Kraft'],
+            'N': ['Zyklenzähler'],
+            's_vert_1': ['Vertikalweg 1'],
+            's_vert_2': ['Vertikalweg 2'],
+            's_vert_3': ['Vertikalweg 3'],
+            's_piston': ['Kolbenweg'],
+        }
+
+    def _process_data(self):
+
+        meta, data = read_geosys(self.data, '015', metadata_ids=['001','003', '005'])
+
+        #define in class
+        self.data = data.reset_index()
+        self.metadata.update(meta)
+
+        # log infos
+        self._logger.debug(f'metadata: {self.metadata}')
+        self._logger.debug(f'data: {self.data.head()}')
+        
+        print(data.head())