add new example

examples/sogno_demo/Sogno_cigre_se.py

+#!/usr/bin/env python
+
+################
+# Dependancies:
+# paho, mqtt library for python (sudo pip install paho-mqtt)
+import paho.mqtt.client as mqtt
+import json
+import numpy as np
+import time
+import logging
+
+logging.basicConfig(filename='recv_client.log', level=logging.INFO)
+
+import sys
+#add to path the path to the state-estimation commit 23facb6b60d3f7993eb05ba7703f137b588343ec
+sys.path.append("../../../state-estimation-client/state-estimation/acs/state_estimation")
+from network import System, Ymatrix_calc
+from nv_state_estimator_cim import DsseCall
+from nv_state_estimator_cim import Complex_to_all
+from nv_state_estimator_cim import Real_to_all
+from measurement_generator_online import Zdata_structure_creation, Zdatameas_creation_fromPF
+import cimpy
+
+Sb = 25e6
+Vb1 = 110e3
+Vb2 = 20e3
+Zb = (Vb2**2)/Sb
+
+def on_message(mqttc, userdata, msg):
+    global end
+    #stop simulation when message "end" is received
+    if msg.payload=="end":
+        end = True
+    else:   
+        message = json.loads(msg.payload)[0]
+        values = message["data"]
+        values_se = np.zeros(len(values))
+        Vmag_est = np.zeros(len(system.nodes))
+        Vphase_est = np.zeros(len(system.nodes))
+        Vmag_true = np.zeros(len(system.nodes))
+        for elem in range(len(system.nodes)):
+            #if int(system.nodes[elem].uuid[1:]) == 0:
+            if system.nodes[elem].uuid == "N0":
+                values_se[elem*2] = values[int(system.nodes[elem].uuid[1:])*2]/Vb1
+            else:
+                values_se[elem * 2] = values[int(system.nodes[elem].uuid[1:])*2]/Vb2
+            values_se[elem*2+1] = values[int(system.nodes[elem].uuid[1:])*2+1]
+    
+        if int(message["sequence"]) > 0:
+            # To delete later
+            nodes_num = len(system.nodes)
+            branches_num = len(system.branches)
+            Vt = np.zeros(len(system.nodes), dtype=np.complex)
+            #index = np.linspace(0, len(values) - 2, int(len(values) / 2)).astype(int)     
+            for elem in range(len(system.nodes)):
+                index = 2*system.nodes[elem].index
+                #Vt[elem] = values_se[index[elem]]*(np.cos(values_se[index[elem]+1]) + 1j*np.sin(values_se[index[elem]+1]))
+                Vt[elem] = values_se[index]*(np.cos(values_se[index+1]) + 1j*np.sin(values_se[index+1]))
+            Vtrue = Complex_to_all(Vt)
+            #print(Vtrue.mag)
+
+            """ From here on, all the other quantities of the grid are calculated """
+            Irx = np.zeros(branches_num, dtype=np.complex)
+            for index in range(branches_num):
+                fr = system.branches[index].start_node.index  # branch.start[idx]-1
+                to = system.branches[index].end_node.index  # branch.end[idx]-1
+                Irx[index] = - (Vtrue.complex[fr] - Vtrue.complex[to]) * Ymatr[fr][to]
+            Ir = np.real(Irx)
+            Ix = np.imag(Irx)
+
+            Itrue = Real_to_all(Ir, Ix)
+            Iinj_r = np.zeros(nodes_num)
+            Iinj_x = np.zeros(nodes_num)
+            for k in range(nodes_num):
+                to = []
+                fr = []
+                for m in range(branches_num):
+                    if k == system.branches[m].start_node.index:
+                        fr.append(m)
+                    if k == system.branches[m].end_node.index:
+                        to.append(m)
+
+                Iinj_r[k] = np.sum(Itrue.real[to]) - np.sum(Itrue.real[fr])
+                Iinj_x[k] = np.sum(Itrue.imag[to]) - np.sum(Itrue.imag[fr])
+
+            Iinjtrue = Real_to_all(Iinj_r, Iinj_x)
+            Sinj_rx = np.multiply(Vtrue.complex, np.conj(Iinjtrue.complex))
+            Sinjtrue = Real_to_all(np.real(Sinj_rx), np.imag(Sinj_rx))
+            values_se = np.append(values_se, Sinjtrue.real)
+            values_se = np.append(values_se, Sinjtrue.imag)
+
+            # till here
+
+            if message["sequence"] < 90:
+                zdatameas = Zdatameas_creation_fromPF(meas_config1, zdata1, values_se, message["sequence"])
+            else:
+                zdatameas = Zdatameas_creation_fromPF(meas_config2, zdata2, values_se, message["sequence"])
+
+            Vest, Iest, Iinjest, S1est, S2est, Sinjest = DsseCall(system, zdatameas, Ymatr, Adj)
+            print("state_estimation_res finish!")
+            
+            Vmag_est = []
+            Vphase_est = []
+            for elem in range(len(system.nodes)):
+                Vmag_est.append(Vest.mag[elem])
+                Vphase_est.append(Vest.phase[elem])
+
+            #send results to dummy_simulator
+            payload = {}
+            payload["client"] = "Sogno_cigre_se"
+            payload["sequence"] = message["sequence"]
+            payload["V_est_mag"] = Vmag_est
+            payload["V_est_phase"] = Vphase_est
+
+            mqttc.publish(topic_dummy_simulator, "[" + json.dumps(payload) + "]", 0)
+
+def on_connect(mqttc, userdata, flags, rc):
+    if rc == 0:
+        mqttc.connected_flag = False  # set flag
+        print("connected OK with returned code=", rc)
+    else:
+        print("Bad connection, error code= ", rc)
+
+xml_files = [
+    r"..\..\..\cim-grid-data\CIGRE_MV\CIGRE_MV_no_tapchanger_With_LoadFlow_Results\Rootnet_FULL_NE_06J16h_EQ.xml",
+    r"..\..\..\cim-grid-data\CIGRE_MV\CIGRE_MV_no_tapchanger_With_LoadFlow_Results\Rootnet_FULL_NE_06J16h_SV.xml",
+    r"..\..\..\cim-grid-data\CIGRE_MV\CIGRE_MV_no_tapchanger_With_LoadFlow_Results\Rootnet_FULL_NE_06J16h_TP.xml"]
+
+res = cimpy.cimread(xml_files)
+cimpy.setNodes(res)
+cimpy.setPowerTransformerEnd(res)
+system = System()
+system.load_cim_data(res)
+system.branches[12].y = system.branches[12].y*(Vb1**2)/(Vb2**2)
+system.branches[13].y = system.branches[13].y*(Vb1**2)/(Vb2**2)
+Ymatrix, Adj = Ymatrix_calc(system)
+Ymatr = Ymatrix*Zb
+
+meas_configfile = r"C:\Users\Martin\Desktop\hiwi\git\state-estimation-client\sogno_demo\Measurement_config.json"
+meas_configfile2 = r"C:\Users\Martin\Desktop\hiwi\git\state-estimation-client\sogno_demo\Measurement_config2.json"
+meas_configfile3 = r"C:\Users\Martin\Desktop\hiwi\git\state-estimation-client\sogno_demo\Measurement_config3.json"
+
+with open(meas_configfile) as f1:
+    payload = json.load(f1)
+
+with open(meas_configfile2) as f2:
+    meas_config1 = json.load(f2)
+
+with open(meas_configfile3) as f3:
+    meas_config2 = json.load(f3)
+
+zdata1 = Zdata_structure_creation(meas_config1, system)
+zdata2 = Zdata_structure_creation(meas_config2, system)
+
+#global variable used in callback on_message
+end = False
+
+#topic names
+topic = "dpsim-powerflow"
+topic_dummy_simulator = "dpsim-powerflow_results"
+
+broker_adress = "m16.cloudmqtt.com"
+mqtt.Client.connected_flag=False                        #create flag in class
+mqttc = mqtt.Client("SognoDemo_Client2", True)          #create new instance
+mqttc.username_pw_set("ilgtdaqk", "UbNQQjmcUdqq")
+mqttc.on_connect = on_connect                           #attach function to callback
+mqttc.on_message = on_message                           #attach function to callback
+mqttc.connect(broker_adress, 14543)	                    #connect to broker
+mqttc.loop_start()                                      #start loop to process callback
+time.sleep(4)                                           #wait for connection setup to complete
+mqttc.subscribe(topic)
+
+#send notification message notifying that the connection was successful
+mqttc.publish(topic_dummy_simulator, json.dumps("OK"))
+
+while not end:
+    time.sleep(0.01)
+
+mqttc.loop_stop()       #Stop loop 
+mqttc.disconnect()      #disconnect
\ No newline at end of file

examples/sogno_demo/Sogno_cigre_se_cim.py

+import time
+import logging
+import json
+import paho.mqtt.client as mqtt
+import numpy as np
+import cimpy
+
+import sys
+sys.path.append("../../acs/state_estimation")
+from network import System
+from results import Results
+from measurement import Measurents_set
+from nv_state_estimator_cim import DsseCall
+
+logging.basicConfig(filename='recv_client.log', level=logging.INFO)
+
+# The callback for when the client receives a CONNACK response from the server.
+def on_connect(client, userdata, flags, rc):
+    if rc==0:
+        client.connected_flag=True      #set flag
+        print("connected OK with returned code=", rc)
+    else:
+        print("Bad connection with returned code=", rc)
+
+# The callback for when a PUBLISH message is received from the server.
+def on_message(client, userdata, msg):
+    global end
+    #stop simulation when message "end" is received
+    if json.loads(msg.payload) == "end":
+        end = True
+    else:
+        message = json.loads(msg.payload)[0]
+        sequence = message['sequence']
+        data = message['data']
+
+        if sequence>0:
+            #store the recived data in powerflow_results
+            for node in powerflow_results.nodes:
+                magnitude = 0.0
+                phase = 0.0
+                uuid = node.topology_node.uuid
+                for elem_idx, elem_data in enumerate(mapping):
+                    if elem_data[0] == uuid:                 #elem_data[0] == uuid
+                        if elem_data[1] == "mag":            #elem_data[1] = "mag" or "phase"
+                            magnitude = data[elem_idx]
+                            #convert to per unit system
+                            if uuid=="N0":
+                                magnitude /= Vb1
+                            else:
+                                magnitude /= Vb2 
+                        elif elem_data[1] == "phase":
+                            phase = data[elem_idx]
+                node.voltage = magnitude*(np.cos(phase) + 1j*np.sin(phase))
+
+            #calculate quantities I, Iinj, S and Sinj
+            powerflow_results.calculate_all()
+
+            #read measurements from file
+            measurements_set = Measurents_set()
+            if sequence<90:
+                measurements_set.read_measurements_from_file(powerflow_results, meas_configfile1)
+            else:
+                measurements_set.read_measurements_from_file(powerflow_results, meas_configfile2)
+           
+            #calculate the measured values (affected by uncertainty)
+            measurements_set.meas_creation(sequence)
+            #print(measurements_set.getmVal_test())
+
+            #Performs state estimation
+            state_estimation_res = DsseCall(system, measurements_set)
+            print("state_estimation_res finish!")
+
+            #send results to dummy_simulator
+            payload = {}
+            payload["client"] = "Sogno_cigre_se_cim"
+            payload["sequence"] = sequence
+            payload["V_est_mag"] = np.absolute(state_estimation_res.get_voltages()).tolist()
+            payload["V_est_phase"] = np.angle(state_estimation_res.get_voltages()).tolist()
+
+            #print(payload)
+            mqttc.publish(topic_dummy_simulator, "[" + json.dumps(payload) + "]", 0)
+
+#grid files
+xml_files = [
+    r"..\..\..\cim-grid-data\CIGRE_MV\CIGRE_MV_no_tapchanger_With_LoadFlow_Results\Rootnet_FULL_NE_06J16h_EQ.xml",
+    r"..\..\..\cim-grid-data\CIGRE_MV\CIGRE_MV_no_tapchanger_With_LoadFlow_Results\Rootnet_FULL_NE_06J16h_SV.xml",
+    r"..\..\..\cim-grid-data\CIGRE_MV\CIGRE_MV_no_tapchanger_With_LoadFlow_Results\Rootnet_FULL_NE_06J16h_TP.xml"]
+
+#measurements files
+meas_configfile1 = r"..\..\..\state-estimation-client\sogno_demo\Measurement_config2.json"
+meas_configfile2 = r"..\..\..\state-estimation-client\sogno_demo\Measurement_config3.json"  
+
+#per unit basis values
+Sb = 25e6
+Vb1 = 110e3
+Vb2 = 20e3
+Zb = (Vb2**2)/Sb
+
+#load grid
+res = cimpy.cimread(xml_files)
+cimpy.setNodes(res)
+cimpy.setPowerTransformerEnd(res)
+system = System()
+system.load_cim_data(res)
+
+#Ymatrix -> per unit
+system.branches[12].y *= (Vb1**2)/(Vb2**2)
+system.branches[13].y *= (Vb1**2)/(Vb2**2)
+system.Ymatrix_calc()
+system.Ymatrix *= Zb
+
+#create a results object to store the received data
+powerflow_results = Results(system)
+
+#read mapping file and split each line of mapping_file by point e.g.: N0.V.mag -> ["V0", "mag"]
+mapping_file = r"..\..\..\state-estimation-client\sogno_demo\villas_sent_data.conf"
+mapping = []
+with open(mapping_file) as mfile:
+    for line in mfile:
+        mapping.append(line.strip('\n'))
+for num, elem in enumerate(mapping):
+    uuid, V, type = elem.split(".")
+    mapping[num] = [None]*2
+    mapping[num][0] = uuid
+    mapping[num][1] = type
+
+#global variable used in callback on_message
+end = False
+
+##topic names
+topic = "dpsim-powerflow"                               
+topic_dummy_simulator = "dpsim-powerflow_results"          
+
+broker_adress = "m16.cloudmqtt.com"
+mqtt.Client.connected_flag=False                        #create flag in class
+mqttc = mqtt.Client("SognoDemo_Client1", True)          #create new instance
+mqttc.username_pw_set("ilgtdaqk", "UbNQQjmcUdqq")
+mqttc.on_connect = on_connect                           #attach function to callback
+mqttc.on_message = on_message                           #attach function to callback
+mqttc.connect(broker_adress, 14543)                     #connect to broker
+mqttc.loop_start()                                      #start loop to process callback
+time.sleep(4)                                           #wait for connection setup to complete
+mqttc.subscribe(topic)
+
+while not mqttc.connected_flag:     #wait in loop
+    print("In wait loop")
+    time.sleep(1)
+
+#send notification message notifying that the connection was successful
+mqttc.publish(topic_dummy_simulator, json.dumps("OK"))
+
+while not end:
+    time.sleep(0.01)
+
+mqttc.loop_stop()       #Stop loop 
+mqttc.disconnect()      #disconnect
\ No newline at end of file

examples/sogno_demo/dummy_simulator.py

+import json
+import paho.mqtt.client as mqtt
+import time
+import numpy as np
+
+# The callback for when the client receives a CONNACK response from the server.
+def on_connect(client, userdata, flags, rc):
+    if rc==0:
+        client.connected_flag=True      #set flag
+        print("connected OK with returned code=",rc)
+    else:
+        print("Bad connection with returned code=",rc)
+
+# The callback for when a PUBLISH message is received from the server.
+def on_message(client, userdata, msg):
+    global connected_clients, sequence, V_est_sogno_cigre_se, V_est_sogno_cigre_se_cim, end
+    if connected_clients<2:
+        #check that the two clients are connected
+        if json.loads(msg.payload) == "OK":
+            connected_clients +=1
+            print("client" + str(connected_clients) + " connected")
+            if connected_clients==2:
+                mqttc.publish(topic_clients, data[sequence])
+    else:
+        #receive and compare the results of the state estimation
+        message = json.loads(msg.payload)[0]
+        from_client = message['client']
+        V_est_mag = np.array(message['V_est_mag'])
+        V_est_phase = np.array(message['V_est_phase'])
+        if from_client == "Sogno_cigre_se_cim":
+            V_est_sogno_cigre_se_cim = V_est_mag + 1j*V_est_phase
+        elif from_client == "Sogno_cigre_se":
+            V_est_sogno_cigre_se = V_est_mag + 1j*V_est_phase
+        if V_est_sogno_cigre_se_cim is not None and V_est_sogno_cigre_se is not None:
+            if (np.array_equal(np.around(V_est_sogno_cigre_se_cim,3), np.around(V_est_sogno_cigre_se,3))):
+                print("V_est_sogno_cigre_se==V_est_sogno_cigre_se_cim (sequence = " + str(sequence) +")?: True ")
+            else:
+                print("V_est_sogno_cigre_se==V_est_sogno_cigre_se_cim (sequence = " + str(sequence) +")?: False ")
+                print(V_est_sogno_cigre_se_cim - V_est_sogno_cigre_se)
+            sequence += 1
+            if sequence == len_data+1:
+                end = True
+                mqttc.publish(topic_clients, json.dumps("end"))
+            else:
+                mqttc.publish(topic_clients, data[sequence])
+                V_est_sogno_cigre_se = None
+                V_est_sogno_cigre_se_cim = None
+
+#global variables used in callback on_message
+end = False
+connected_clients = 0 
+sequence = 1
+V_est_sogno_cigre_se = None
+V_est_sogno_cigre_se_cim = None
+len_data = 300          #number of sequences to use
+
+#topic names
+topic = "dpsim-powerflow_results"
+topic_clients = "dpsim-powerflow"
+
+broker_adress = "m16.cloudmqtt.com"
+mqttc = mqtt.Client("SognoDemo", True)          #create new instance
+mqttc.username_pw_set("ilgtdaqk", "UbNQQjmcUdqq")
+mqttc.on_connect = on_connect                   #attach function to callback
+mqttc.on_message = on_message                   #attach function to callback
+mqttc.connect(broker_adress, 14543)             #connect to broker
+mqttc.loop_start()                              #start loop to process callback
+time.sleep(1)                                   #wait for connection setup to complete
+mqttc.subscribe(topic)
+
+data_file = r"C:\Users\Martin\Desktop\hiwi\git\state-estimation-client\sogno_demo\dpsim_powerflow_record_cigre.txt"
+data = []
+with open(data_file) as json_file:
+    for line in json_file:
+        data.append(line)
+
+while not end:
+    time.sleep(0.01)
+
+mqttc.loop_stop()       #Stop loop 
+mqttc.disconnect()      #disconnect
\ No newline at end of file

examples/sogno_demo/init.bat

+@echo off
+
+set PYTHON_DIR=C:\Python36
+set SOGNO_CIGRE_SE=Sogno_cigre_se.py
+set SOGNO_CIGRE_SE_CIM=Sogno_cigre_se_cim.py
+set DUMMY_SIMULATOR=dummy_simulator.py
+
+
+start CMD /K "title Dummy_Simulator & python %DUMMY_SIMULATOR%"
+TIMEOUT /T 3
+start CMD /K "title Sogno_cigre_se_cim & python %SOGNO_CIGRE_SE_CIM%"
+start CMD /K "title Sogno_cigre_se & python %SOGNO_CIGRE_SE%"