Commit 440596b1 authored by martin.moraga's avatar martin.moraga
Browse files

fixed some errors

parent 7002634e
......@@ -41,10 +41,13 @@ class Measurement():
self.element_type = element_type
self.meas_type = meas_type
self.meas_value = meas_value
self.std_dev = unc/100
"""
if meas_type not in [MeasType.Ipmu_phase, MeasType.Vpmu_phase]:
self.std_dev = meas_value*unc/100
elif meas_type in [MeasType.Ipmu_phase, MeasType.Vpmu_phase]:
self.std_dev = unc/100
"""
self.mval = 0.0 #measured values (affected by uncertainty)
class Measurents_set():
......@@ -66,7 +69,6 @@ class Measurents_set():
"""
with open(file_name) as json_file:
data = json.load(json_file)
for key, value in data['Measurement'].items():
if key=="Vmag":
unc = float(value['unc'])
......@@ -122,8 +124,12 @@ class Measurents_set():
for index in value['idx']:
element = powerflow_results.nodes[index-1].topology_node
meas_value_mag = np.abs(powerflow_results.nodes[index-1].voltage)
meas_value_phase = np.angle(powerflow_results.nodes[index-1].voltage)
#meas_value_phase = np.angle(powerflow_results.nodes[index-1].voltage)
self.create_measurement(element, ElemType.Node, MeasType.Vpmu_mag, meas_value_mag, unc_mag)
#self.create_measurement(element, ElemType.Node, MeasType.Vpmu_phase, meas_value_phase, unc_phase)
for index in value['idx']:
element = powerflow_results.nodes[index-1].topology_node
meas_value_phase = np.angle(powerflow_results.nodes[index-1].voltage)
self.create_measurement(element, ElemType.Node, MeasType.Vpmu_phase, meas_value_phase, unc_phase)
elif key=="Ipmu":
unc_mag = float(value['unc_mag'])
......@@ -151,7 +157,12 @@ class Measurents_set():
err_pu = np.random.uniform(-1,1,len(self.measurements))
for index, measurement in enumerate(self.measurements):
measurement.mval = measurement.meas_value + measurement.std_dev*err_pu[index]
if measurement.meas_type not in [MeasType.Ipmu_phase, MeasType.Vpmu_phase]:
zdev = measurement.meas_value*measurement.std_dev
elif measurement.meas_type in [MeasType.Ipmu_phase, MeasType.Vpmu_phase]:
zdev = measurement.std_dev
measurement.mval = measurement.meas_value + zdev*err_pu[index]
def meas_creation_test(self, err_pu):
"""
......@@ -206,7 +217,7 @@ class Measurents_set():
#the weight is small and can bring instability during matrix inversion, so we "cut" everything below 10^-6
if measurement.std_dev<10**(-6):
measurement.std_dev = 10**(-6)
weights[index] = measurement.std_dev**(-2)
weights[index] = (measurement.std_dev/3)**(-2)
return weights
......@@ -262,4 +273,14 @@ class Measurents_set():
for index, measurement in enumerate(self.measurements):
std_dev[index] = measurement.std_dev
return std_dev
\ No newline at end of file
return std_dev
def getmVal_test(self):
"""
returns an array with all measured values (affected by uncertainty)
"""
mVal = np.zeros(len(self.measurements))
for index, measurement in enumerate(self.measurements):
mVal[index] = measurement.mval
return mVal
\ No newline at end of file
......@@ -18,10 +18,10 @@ def DsseCall(system, measurements):
Vmag_meas = 0
Vpmu_meas = 0
for elem in measurements.measurements:
if elem.meas_type=="V_mag":
Vmag_meas = Vmag_meas + 1
elif elem.meas_type=="Vpmu_mag":
Vpmu_meas = Vpmu_meas + 1
if elem.meas_type==MeasType.V_mag:
Vmag_meas += 1
elif elem.meas_type==MeasType.Vpmu_mag:
Vpmu_meas += 1
trad_code = 1 if Vmag_meas >0 else 0
PMU_code = 2 if Vpmu_meas >0 else 0
......@@ -165,11 +165,10 @@ def DssePmu(nodes_num, measurements, Gmatrix, Bmatrix, Adj):
@param Adj
return: np.array V - estimated voltages
"""
#calculate weights matrix (obtained as stdandard_deviations^-2)
weights = measurements.getWeightsMatrix()
W = np.diag(weights)
#Jacobian for Power Injection Measurements
H2, H3 = calculateJacobiMatrixSinj(measurements, nodes_num, Gmatrix, Bmatrix, Adj, type=1)
......@@ -191,7 +190,7 @@ def DssePmu(nodes_num, measurements, Gmatrix, Bmatrix, Adj):
WH = np.inner(W,H.transpose())
G = np.inner(H.transpose(),WH.transpose())
Ginv = np.linalg.inv(G)
#get array which contains the index of measurements type MeasType.Sinj_real, MeasType.Sinj_imag in the array measurements.measurements
pidx = measurements.getIndexOfMeasurements(type=MeasType.Sinj_real)
qidx = measurements.getIndexOfMeasurements(type=MeasType.Sinj_imag)
......@@ -204,12 +203,12 @@ def DssePmu(nodes_num, measurements, Gmatrix, Bmatrix, Adj):
#get an array with all measured values (affected by uncertainty)
z = measurements.getmVal()
V = np.ones(nodes_num) + 1j*np.zeros(nodes_num)
State = np.concatenate((np.ones(nodes_num), np.zeros(nodes_num)), axis=0)
epsilon = 5
num_iter = 0
# Iteration of Netwon Rapson method: needed to solve non-linear system of equation
while epsilon>10**(-6):
......@@ -217,7 +216,7 @@ def DssePmu(nodes_num, measurements, Gmatrix, Bmatrix, Adj):
# in every iteration the input power measurements are converted into currents by dividing by the voltage estimated at the previous iteration
z = convertSinjMeasIntoCurrents(measurements, V, z, pidx, qidx)
z = convertSbranchMeasIntoCurrents(measurements, V, z, p1br, q1br, p2br, q2br)
""" WLS computation """
y = np.inner(H,State)
res = np.subtract(z,y)
......@@ -231,6 +230,7 @@ def DssePmu(nodes_num, measurements, Gmatrix, Bmatrix, Adj):
V.real = State[:nodes_num]
V.imag = State[nodes_num:]
#print(V.real)
num_iter = num_iter+1
return V
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment