Merge branch 'master' of git.rwth-aachen.de:acs/core/simulation/data-processing

dataprocessing/calc.py

-import matplotlib.pyplot as plt
-import numpy as np
-from .timeseries import *
-
-def diff(name, ts1, ts2):
-    """ Calculate difference.
-    Assumes the same time steps for both timeseries.
-    """
-    ts_diff = TimeSeries(name, ts1.time, (ts1.values - ts2.values))
-    return ts_diff
-
-def scale_ts(name, ts, factor):
-    """ Scale timeseries.
-    Assumes the same time steps for both timeseries.
-    """
-    ts_scaled = TimeSeries(name, ts.time, ts.values * factor)
-    return ts_scaled
-
-def complex_abs(name, real, imag):
-    """ Calculate absolute value of complex variable.
-    Assumes the same time steps for both timeseries.
-    """
-    ts_abs = TimeSeries(name, real.time, np.sqrt(real.values ** 2 + imag.values ** 2))
-    return ts_abs
-

dataprocessing/plotdpsim.py

@@ -6,8 +6,8 @@ import matplotlib.pyplot as plt
 matplotlib.rcParams.update({'font.size': 8})
 
 def plot_dpsim_abs_diff(filename1, label1, node1, filename2, label2, node2):
-    ts_dpsim1 = read_time_series_DPsim(filename1)
-    ts_dpsim2 = read_time_series_DPsim(filename2)
+    ts_dpsim1 = read_timeseries_DPsim(filename1)
+    ts_dpsim2 = read_timeseries_DPsim(filename2)
 
     ts_dpsim1_length = len(ts_dpsim1)
     im_offset1 = int(ts_dpsim1_length / 2)
@@ -44,8 +44,8 @@ def plot_dpsim_abs_diff(filename1, label1, node1, filename2, label2, node2):
     plt.show()
 
 def plot_dpsim_abs(filename1, label1, node1, filename2, label2, node2):
-        ts_dpsim1 = read_time_series_DPsim(filename1)
-        ts_dpsim2 = read_time_series_DPsim(filename2)
+        ts_dpsim1 = read_timeseries_DPsim(filename1)
+        ts_dpsim2 = read_timeseries_DPsim(filename2)
 
         ts_dpsim1_length = len(ts_dpsim1)
         im_offset1 = int(ts_dpsim1_length / 2)
@@ -79,8 +79,47 @@ def plot_dpsim_abs(filename1, label1, node1, filename2, label2, node2):
         plt.show()
 
 
+def plot_dpsim_emt_abs(filenameDP, nodeDP, filenameEMT, nodeEMT):
+    ts_dpsimDP = read_timeseries_DPsim(filenameDP)
+    ts_dpsimEMT = read_timeseries_DPsim(filenameEMT)
+
+    ts_dpsimDP_length = len(ts_dpsimDP)
+    im_offsetDP = int(ts_dpsimDP_length / 2)
+    if im_offsetDP <= nodeDP or nodeDP < 0:
+        print('Node DP not available')
+        exit()
+
+    ts_dpsimEMT_length = len(ts_dpsimEMT)
+    if ts_dpsimEMT_length <= nodeEMT or nodeEMT < 0:
+        print('Node EMT not available')
+        exit()
+
+    ts_absDP = complex_abs('node ' + str(nodeDP) + 'abs', ts_dpsimDP[nodeDP], ts_dpsimDP[nodeDP + im_offsetDP])
+    ts_absDP = scale_ts(ts_absDP.name, ts_absDP, 0.001)
+    ts_absDP.label = 'DP abs'
+
+    ts_shiftDP = dyn_phasor_shift_to_emt('node ' + str(nodeDP) + 'shift', ts_dpsimDP[nodeDP], ts_dpsimDP[nodeDP + im_offsetDP], 50)
+    ts_shiftDP = scale_ts(ts_shiftDP.name, ts_shiftDP, 0.001)
+    ts_shiftDP.label = 'DP shift'
+
+    ts_EMT = TimeSeries('node ' + str(nodeEMT), ts_dpsimEMT[nodeEMT].time, ts_dpsimEMT[nodeEMT].values)
+    ts_EMT = scale_ts(ts_EMT.name, ts_EMT, 0.001)
+    ts_EMT.label = 'EMT'
+
+    figure_id = 1
+    # plt.figure(figure_id)
+    plt.figure(figure_id, figsize=(12 / 2.54, 6 / 2.54), facecolor='w', edgecolor='k')
+    plot_timeseries(figure_id, ts_EMT)
+    plot_timeseries(figure_id, ts_absDP)
+    plot_timeseries(figure_id, ts_shiftDP)
+    plt.xlabel('Time [s]')
+    plt.ylabel('Voltage [kV]')
+    plt.grid(True)
+    plt.tight_layout()
+    plt.show()
+
 def plot_dpsim_abs_single(filename, node):
-    ts_dpsim = read_time_series_DPsim(filename)
+    ts_dpsim = read_timeseries_DPsim(filename)
 
     ts_dpsim_length = len(ts_dpsim)
     print('DPsim results file length:')

dataprocessing/plottools.py

@@ -16,7 +16,7 @@ def plot_timeseries(figure_id, timeseries, plt_linestyle='-', plt_linewidth=2, p
         else:
             plt.plot(timeseries.time, timeseries.values, linestyle=plt_linestyle, label=timeseries.label, linewidth=plt_linewidth)
         plt.gca().autoscale(axis='x', tight=True)
-        plt.legend()
+        plt.legend(loc='lower right')
     else:
         for ts in timeseries:
             plt.subplot(len(timeseries), 1, timeseries.index(ts) + 1)

dataprocessing/readtools.py

@@ -35,16 +35,38 @@ def read_timeseries_PLECS(filename, timeseries_names=None):
             timeseries_list.append(TimeSeries(name, pd_df['Time'].values, pd_df[name].values))
     return timeseries_list
 
-def read_timeseries_DPsim(filename, timeseries_names=None):
-    pd_df = pd.read_csv(filename)
+def read_timeseries_dpsim_real(filename, header=None, timeseries_names=None):
+    """Reads real time series data from DPsim log file which may have a header.
+    Timeseries names are assigned according to the header names if available.
+    :param filename: name of the csv file that has the data
+    :param header: specifies if the log file has a header
+    :param timeseries_names: column names which should be read
+    :return: list of Timeseries objects
+    """
     timeseries_list = []
 
+    if header is True:
+        pd_df = pd.read_csv(filename)
+    else:
+        pd_df = pd.read_csv(filename, header=None)
+
     if timeseries_names is None:
         # No trajectory names specified, thus read in all
-        timeseries_names = list(pd_df.columns.values)
-        timeseries_names.remove('Time')
-        for name in timeseries_names:
-            timeseries_list.append(TimeSeries(name, pd_df['Time'].values, pd_df[name].values))
+        column_names = list(pd_df.columns.values)
+        # Remove timestamps column name and store separately
+        column_names.remove(0)
+        timestamps = pd_df.iloc[:,0]
+
+        if header is True:
+            for name in column_names:
+                timeseries_list.append(TimeSeries(name, timestamps, pd_df[name].values))
+        else:
+            node_number = int(len(column_names))
+            node_index = 1
+            for column in column_names:
+                ts_name = 'node ' + str(node_index)
+                timeseries_list.append(TimeSeries(ts_name, timestamps, pd_df.iloc[:, column]))
+                node_index = node_index + 1
     else:
         # Read in specified time series
         print('no column names specified yet')
@@ -55,23 +77,68 @@ def read_timeseries_DPsim(filename, timeseries_names=None):
         print(result.name)
     return timeseries_list
 
-def read_timeseries_DPsim_node_values(filename, timeseries_names=None):
+def read_timeseries_dpsim_cmpl(filename, timeseries_names=None):
+    """Reads complex time series data from DPsim log file. Real and
+    imaginary part are stored in one complex variable.
+    :param filename: name of the csv file that has the data
+    :param timeseries_names: column name which should be read
+    :return: list of Timeseries objects
+    """
     pd_df = pd.read_csv(filename, header=None)
     timeseries_list = []
 
     if timeseries_names is None:
         # No trajectory names specified, thus read in all
         column_names = list(pd_df.columns.values)
+        # Remove timestamps column name and store separately
         column_names.remove(0)
+        timestamps = pd_df.iloc[:,0]
+        # Calculate number of network nodes since array is [real, imag]
+        node_number = int(len(column_names) / 2)
         node_index = 1
+        for column in column_names:
+            if node_index <= node_number:
+                ts_name = 'node '+ str(node_index)
+                timeseries_list.append(TimeSeries(ts_name, timestamps, np.vectorize(complex)(pd_df.iloc[:,column],pd_df.iloc[:,column + node_number])))
+            else:
+                break
+            node_index = node_index + 1
+    else:
+        # Read in specified time series
+        print('cannot read specified columns yet')
+
+    print('DPsim results file length:')
+    print(len(timeseries_list))
+    for result in timeseries_list:
+        print(result.name)
+    return timeseries_list
+
+def read_timeseries_dpsim_cmpl_separate(filename, timeseries_names=None):
+    """Deprecated - Reads complex time series data from DPsim log file. Real and
+    imaginary part are stored separately.
+    :param filename: name of the csv file that has the data
+    :param timeseries_names: column name which should be read
+    :return: list of Timeseries objects
+    """
+    pd_df = pd.read_csv(filename, header=None)
+    timeseries_list = []
+
+    if timeseries_names is None:
+        # No trajectory names specified, thus read in all
+        column_names = list(pd_df.columns.values)
+        # Remove timestamps column name and store separately
+        column_names.remove(0)
+        timestamps = pd_df.iloc[:, 0]
+        # Calculate number of network nodes since array is [real, imag]
         node_number = int(len(column_names) / 2)
+        node_index = 1
         for column in column_names:
             if node_index <= node_number:
-                node_name = node_index
-                timeseries_list.append(TimeSeries('node '+ str(node_name) +' Re', pd_df.iloc[:,0], pd_df.iloc[:,column]))
+                node_name = 'node '+ str(node_index) +' Re'
+                timeseries_list.append(TimeSeries(node_name, timestamps, pd_df.iloc[:,column]))
             else:
-                node_name = node_index - node_number
-                timeseries_list.append(TimeSeries('node '+ str(node_name) +' Im', pd_df.iloc[:,0], pd_df.iloc[:,column]))
+                node_name = 'node '+ str(node_index - node_number) +' Im'
+                timeseries_list.append(TimeSeries(node_name, timestamps, pd_df.iloc[:,column]))
 
             node_index = node_index + 1
     else:

dataprocessing/timeseries.py

 import numpy as np
 
 class TimeSeries:
+    """Stores data from different simulation sources.
+    A TimeSeries object always consists of timestamps and datapoints.
+    """
     def __init__(self, name, time, values, label=""):
         self.time = np.array(time)
         self.values = np.array(values)
         self.name = name
-        self.label = name
\ No newline at end of file
+        self.label = name
+
+    @staticmethod
+    def diff(name, ts1, ts2):
+        """Returns difference between values of two Timeseries objects.
+        Assumes the same time steps for both timeseries.
+        """
+        ts_diff = TimeSeries(name, ts1.time, (ts1.values - ts2.values))
+        return ts_diff
+
+
+    def scale_ts(self, name, factor):
+        """Returns scaled timeseries.
+        Assumes the same time steps for both timeseries.
+        """
+        ts_scaled = TimeSeries(name, self.time, self.values * factor)
+        return ts_scaled
+
+    @staticmethod
+    def complex_abs_dep(name, ts_real, ts_imag):
+        """ Calculate absolute value of complex variable.
+        Assumes the same time steps for both timeseries.
+        """
+        ts_abs = TimeSeries(name, ts_real.time, np.sqrt(ts_real.values ** 2 + ts_imag.values ** 2))
+        return ts_abs
+
+    @staticmethod
+    def complex_abs(name, ts_real, ts_imag):
+        """ Calculate absolute value of complex variable.
+        Assumes the same time steps for both timeseries.
+        """
+        ts_complex = np.vectorize(complex)(ts_real.values, ts_imag.values)
+        ts_abs = TimeSeries(name, ts_real.time, ts_complex.abs())
+        return ts_abs
+
+    def abs(self, name):
+        """ Calculate absolute value of complex variable.
+        Assumes the same time steps for both timeseries.
+        """
+        ts_abs = TimeSeries(name, self.time, self.values.abs())
+        return ts_abs
+
+    def complex_phase(name, ts_real, ts_imag):
+        """ Calculate absolute value of complex variable.
+        Assumes the same time steps for both timeseries.
+        """
+        ts_complex = np.vectorize(complex)(ts_real.values, ts_imag.values)
+        ts_abs = TimeSeries(name, ts_real.time, ts_complex.phase())
+        return ts_abs
+
+    @staticmethod
+    def dyn_phasor_shift_to_emt(name, real, imag, freq):
+        """ Shift dynamic phasor values to EMT by frequency freq.
+            Assumes the same time steps for both timeseries.
+        """
+        ts_shift = TimeSeries(name, real.time, real.values*np.cos(2*np.pi*freq*real.time) - imag.values*np.sin(2*np.pi*freq*real.time))
+        return ts_shift
+
+    @staticmethod
+    def check_node_number_comp(ts_comp, node):
+        """
+        Check if node number is available in complex time series.
+        :param ts_comp: complex time series
+        :param node: node number to be checked
+        :return: true if node number is available, false if out of range
+        """
+        ts_comp_length = len(ts_comp)
+        im_offset = int(ts_comp_length / 2)
+        if im_offset <= node or node < 0:
+            print('Complex node not available')
+            return false
+        else:
+            return true
+
+    @staticmethod
+    def check_node_number(ts, node):
+        """
+        Check if node number is available in time series.
+        :param ts: time series
+        :param node: node number to be checked
+        :return: true if node number is available, false if out of range
+        """
+        ts_length = len(ts)
+        if ts_length <= node or node < 0:
+            print('Node not available')
+            return false
+        else:
+            return true
\ No newline at end of file