Add basic uncert values

This commit is contained in:
Cal Wing 2024-10-22 15:46:08 +10:00
parent b86694e603
commit 578378f23e
9 changed files with 119 additions and 20 deletions

View File

@ -49,6 +49,7 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 3.05 # mm
gauge-c2c: 4 #mm
dist-uncert: 2 #mm
gauges:
- 1
@ -60,6 +61,7 @@ probe-info:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

View File

@ -49,6 +49,7 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 0.8 # mm
gauge-c2c: 1.8 #mm
dist-uncert: 0.025 #mm
gauges:
- 1
@ -61,6 +62,7 @@ probe-info:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

View File

@ -52,6 +52,7 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 0.8 # mm
gauge-c2c: 1.8 #mm
dist-uncert: 0.025 #mm
gauges:
- 1
@ -64,6 +65,7 @@ probe-info:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

View File

@ -54,6 +54,7 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 0.8 # mm
gauge-c2c: 1.8 #mm
dist-uncert: 0.025 #mm
gauges:
- 1
@ -66,6 +67,7 @@ probe-info:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

View File

@ -52,6 +52,7 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 0.8 # mm
gauge-c2c: 1.8 #mm
dist-uncert: 0.025 #mm
gauges:
- 1
@ -64,6 +65,7 @@ probe-info:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

View File

@ -55,11 +55,13 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 0.8 # mm
gauge-c2c: 1.8 #mm
dist-uncert: 0.025 #mm
data-records:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

View File

@ -58,11 +58,13 @@ probe-info:
c2c: 5.6 # mm
gauge-diam: 0.8 # mm
gauge-c2c: 1.8 #mm
dist-uncert: 0.025 #mm
data-records:
- type: "scope"
config: "eProbe-Scope.txt"
data: "eProbe-Scope.csv"
time-uncert: 0 # s
trigger: # Redundant?
type: "channel"

91
main.py
View File

@ -39,6 +39,21 @@ data_to_load = [
"x2s5832"
]
# ==== Uncerts ====
# Taken from DOI: 10.1007/s00193-017-0763-3 (Implementation of a state-to-state analytical framework for the calculation of expansion tube flow properties)
UNCERTS = TUNNEL_INFO["uncertainties"]
def deltaX(delta_x_1: float, delta_x_2: float):
return np.sqrt(np.pow(delta_x_1, 2) + np.pow(delta_x_2, 2))
def deltaT(delta_t_1: float, delta_t_2: float, delta_t_sr: float):
return np.sqrt(np.pow(delta_t_1, 2) + np.pow(delta_t_2, 2) + np.pow(delta_t_sr, 2))
def deltaVs(V: float, dx: float, dt: float, delta_x: tuple[float, float], delta_t: tuple[float, float, float]):
return V * np.sqrt(np.pow(deltaX(*delta_x) / dx, 2) + np.pow(deltaT(*delta_t) / dt, 2))
# ==== Data Loading & Processing ====
def load_data(data_path: str, data={}) -> dict:
data_info_path = data_path + DATA_INFO
@ -114,6 +129,7 @@ def load_data(data_path: str, data={}) -> dict:
"x2": None,
"probes": None, # This may be x2 but may not - ie a scope was used
"trigger_index": None,
"probe_uncert": None, #s
},
"data": {
"x2": {}, # Only pop channels with a voltage scale in ./tunnel-info.yaml
@ -196,6 +212,7 @@ def load_data(data_path: str, data={}) -> dict:
if "scope" in data_locs:
data[x2_shot]["data"]["probes"] = [data[x2_shot]["data"]["scope"][1], data[x2_shot]["data"]["scope"][2]]
data[x2_shot]["time"]["probes"] = data[x2_shot]["time"]["scope"]
data[x2_shot]["time"]["probe_uncert"] = scope_data_info["time-uncert"]
# Find Shock Times
@ -215,7 +232,7 @@ def load_data(data_path: str, data={}) -> dict:
first_value, first_value_uncertainty, _, _ = canny_shock_finder(x2_time_us, refData, sigma=sigma, post_suppression_threshold=post_sup_thresh, plot=False, print_func=None)
shock_point = np.where(x2_time_us >= first_value)[0][0] # [BUG] Seems to give n+1
data[x2_shot]["shock-point"][ref] = shock_point, first_value
data[x2_shot]["shock-point"][ref] = shock_point, first_value, first_value_uncertainty
@ -250,7 +267,7 @@ def load_data(data_path: str, data={}) -> dict:
shock_point = np.where(scope_time >= first_value)[0][0] # [BUG] Seems to give n+1
data[x2_shot]["shock-point"][f"{probe}-g1"] = shock_point, first_value
data[x2_shot]["shock-point"][f"{probe}-g1"] = shock_point, first_value, first_value_uncertainty
if 2 in dataInfo["probe-info"]["gauges"]:
# Do the same for G2
@ -264,7 +281,7 @@ def load_data(data_path: str, data={}) -> dict:
raise ValueError(f"{probe}-g2 not detected")
shock_point = np.where(scope_time >= first_value)[0][0] # [BUG] Seems to give n+1
data[x2_shot]["shock-point"][f"{probe}-g2"] = shock_point, first_value
data[x2_shot]["shock-point"][f"{probe}-g2"] = shock_point, first_value, first_value_uncertainty
@ -276,11 +293,19 @@ def load_data(data_path: str, data={}) -> dict:
if i == 0: continue
p1_time = data[x2_shot]["shock-point"][refProbe][1] / 1e6 # Convert to seconds
p2_time = data[x2_shot]["shock-point"][dataInfo["pcb-refs"][i-1]][1] / 1e6 # Convert to seconds
p2p_dist = (TUNNEL_INFO["distance"][refProbe] - TUNNEL_INFO["distance"][dataInfo["pcb-refs"][i-1]]) / 1000 # convert to m
probe_velocity = p2p_dist / abs(p2_time - p1_time) # m/s
print(f"{refProbe}-{dataInfo["pcb-refs"][i-1]} Measured a shock speed of {probe_velocity:.2f} m/s ({probe_velocity/1000:.2f} km/s)")
data[x2_shot]["shock-speed"][f"{refProbe}-{dataInfo["pcb-refs"][i-1]}"] = probe_velocity, True # Speed, Ref
p2p_dist = (TUNNEL_INFO["distance"][refProbe] - TUNNEL_INFO["distance"][dataInfo["pcb-refs"][i-1]]) / 1000 # convert to m
p2p_time = abs(p2_time - p1_time)
probe_velocity = p2p_dist / p2p_time # m/s
p1_time_uncert = data[x2_shot]["shock-point"][refProbe][2] / 1e6 # Convert to seconds
p2_time_uncert = data[x2_shot]["shock-point"][dataInfo["pcb-refs"][i-1]][2] / 1e6 # Convert to seconds
uncert = deltaVs(probe_velocity, p2p_dist, p2p_time, (UNCERTS["probe-dist"][refProbe], UNCERTS["probe-dist"][dataInfo["pcb-refs"][i-1]]), (p1_time_uncert, p2_time_uncert, UNCERTS["time"]["x2-daq"]))
print(f"{refProbe}-{dataInfo["pcb-refs"][i-1]} Measured a shock speed of {probe_velocity:.2f} +/- {uncert:.2f} m/s ({probe_velocity/1000:.2f} +/- {uncert/1000:.2f} km/s [{uncert/probe_velocity * 100 :.2f}%])")
data[x2_shot]["shock-speed"][f"{refProbe}-{dataInfo["pcb-refs"][i-1]}"] = probe_velocity, uncert, True # Speed, Ref
print()
for probe in dataInfo["probe-info"]["locations"]:
@ -288,11 +313,17 @@ def load_data(data_path: str, data={}) -> dict:
g1_time = data[x2_shot]["shock-point"][f"{probe}-g1"][1] / 1e6 # Convert to seconds
g2_time = data[x2_shot]["shock-point"][f"{probe}-g2"][1] / 1e6 # Convert to seconds
c2c_dist = dataInfo["probe-info"]["c2c"] / 1000 # convert to m
c2c_time = abs(g2_time - g1_time)
probe_velocity = c2c_dist / abs(g2_time - g1_time) # m/s
probe_velocity = c2c_dist / c2c_time # m/s
print(f"{probe} Measured a shock speed of {probe_velocity:.2f} m/s ({probe_velocity/1000:.2f} km/s)")
data[x2_shot]["shock-speed"][probe] = probe_velocity, False # Speed, Ref # m/s
g1_time_uncert = data[x2_shot]["shock-point"][f"{probe}-g1"][2] / 1e6 # Convert to seconds
g2_time_uncert = data[x2_shot]["shock-point"][f"{probe}-g2"][2] / 1e6 # Convert to seconds
uncert = deltaVs(probe_velocity, p2p_dist, p2p_time, (0.05/1000, 0.05/1000), (g1_time_uncert, g2_time_uncert, data[x2_shot]["time"]["probe_uncert"]))
print(f"{probe} Measured a shock speed of {probe_velocity:.2f} +/- {uncert:.2f} m/s ({probe_velocity/1000:.2f} +/- {uncert/1000:.2f} km/s)")
data[x2_shot]["shock-speed"][probe] = probe_velocity, uncert, False # Speed, Ref # m/s
else:
print(f"Unable to calculate probe velocity, only have one probe: {f"{probe}-g2" if f"{probe}-g2" in data[x2_shot]["shock-point"] else f"{probe}-g1"}")
@ -300,23 +331,37 @@ def load_data(data_path: str, data={}) -> dict:
for i in range(len(dataInfo["probe-info"]["locations"]) - 1):
probe_locs = dataInfo["probe-info"]["locations"]
p2p = (TUNNEL_INFO["distance"][probe_locs[1]] - TUNNEL_INFO["distance"][probe_locs[0]]) / 1000 # convert to m
p2p_dist = (TUNNEL_INFO["distance"][probe_locs[1]] - TUNNEL_INFO["distance"][probe_locs[0]]) / 1000 # convert to m
if f"{probe_locs[i]}-g1" in data[x2_shot]["shock-point"] and f"{probe_locs[i+1]}-g1" in data[x2_shot]["shock-point"]:
p1_g1_time = data[x2_shot]["shock-point"][f"{probe_locs[i]}-g1"][1] / 1e6 # Convert to seconds
p2_g1_time = data[x2_shot]["shock-point"][f"{probe_locs[i+1]}-g1"][1] / 1e6 # Convert to seconds
p2p_1 = p2p / abs(p2_g1_time - p1_g1_time) # m/s
p2p_time = abs(p2_g1_time - p1_g1_time)
p2p_1 = p2p_dist / p2p_time # m/s
print(f"{probe_locs[i]}-{probe_locs[i + 1]} - G1 - Measured a shock speed of {p2p_1:.2f} m/s ({p2p_1/1000:.2f} km/s)")
data[x2_shot]["shock-speed"][f"{probe_locs[i]}-{probe_locs[i + 1]}-g1"] = p2p_1, False # Speed, Ref
p1_time_uncert = data[x2_shot]["shock-point"][f"{probe_locs[i]}-g1"][2] / 1e6 # Convert to seconds
p2_time_uncert = data[x2_shot]["shock-point"][f"{probe_locs[i+1]}-g1"][2] / 1e6 # Convert to seconds
uncert = deltaVs(p2p_1, p2p_dist, p2p_time, (UNCERTS["probe-dist"][probe_locs[i]], UNCERTS["probe-dist"][probe_locs[i+1]]), (p1_time_uncert, p2_time_uncert, data[x2_shot]["time"]["probe_uncert"]))
print(f"{probe_locs[i]}-{probe_locs[i + 1]} - G1 - Measured a shock speed of {p2p_1:.2f} +/- {uncert:.2f} m/s ({p2p_1/1000:.2f} +/- {uncert/1000:.2f} [{uncert/p2p_1 * 100:.2f}%] km/s)")
data[x2_shot]["shock-speed"][f"{probe_locs[i]}-{probe_locs[i + 1]}-g1"] = p2p_1, uncert, False # Speed, Ref
if f"{probe_locs[i]}-g2" in data[x2_shot]["shock-point"] and f"{probe_locs[i+1]}-g2" in data[x2_shot]["shock-point"]:
p1_g2_time = data[x2_shot]["shock-point"][f"{probe_locs[i]}-g2"][1] / 1e6 # Convert to seconds
p2_g2_time = data[x2_shot]["shock-point"][f"{probe_locs[i+1]}-g2"][1] / 1e6 # Convert to seconds
p2p_2 = p2p / abs(p2_g2_time - p1_g2_time) # m/s
p2p_time = abs(p2_g2_time - p1_g2_time)
p2p_2 = p2p_dist / p2p_time # m/s
print(f"{probe_locs[i]}-{probe_locs[i + 1]} - G2 - Measured a shock speed of {p2p_2:.2f} m/s ({p2p_2/1000:.2f} km/s)")
data[x2_shot]["shock-speed"][f"{probe_locs[i]}-{probe_locs[i + 1]}-g2"] = p2p_2, False # Speed, Ref
p1_time_uncert = data[x2_shot]["shock-point"][f"{probe_locs[i]}-g2"][2] / 1e6 # Convert to seconds
p2_time_uncert = data[x2_shot]["shock-point"][f"{probe_locs[i+1]}-g2"][2] / 1e6 # Convert to seconds
uncert = deltaVs(p2p_2, p2p_dist, p2p_time, (UNCERTS["probe-dist"][probe_locs[i]], UNCERTS["probe-dist"][probe_locs[i+1]]), (p1_time_uncert, p2_time_uncert, data[x2_shot]["time"]["probe_uncert"]))
print(f"{probe_locs[i]}-{probe_locs[i + 1]} - G2 - Measured a shock speed of {p2p_2:.2f} +/- {uncert:.2f} m/s ({p2p_2/1000:.2f} +/- {uncert/1000:.2f} [{uncert/p2p_2 * 100:.2f}%] km/s)")
data[x2_shot]["shock-speed"][f"{probe_locs[i]}-{probe_locs[i + 1]}-g2"] = p2p_2, uncert, False # Speed, Ref
print()
# Return the data & the successfully loaded data keys
@ -344,6 +389,11 @@ def genGraph(gData: dict, showPlot: bool = True, doLimits: bool = True, forcePlo
"plots": []
}
#if forcePlots or not doLimits: graphData["title"] += "\n"
#if forcePlots: graphData["title"] += "(All Data Shown)"
#if not doLimits: graphData["title"] += () + "Full Re"
lims = []
for label in gData["info"]["pcb-refs"]: # + ["trigbox"]:
@ -398,13 +448,14 @@ def genGraph(gData: dict, showPlot: bool = True, doLimits: bool = True, forcePlo
for shock_speed_loc in gData["shock-speed"]:
probeText += "\n"
#probeText += "(Reference) " if gData["shock-speed"][shock_speed_loc][1] else ""
probeText += f"{shock_speed_loc} - {gData["shock-speed"][shock_speed_loc][0]/1000:.2f} km/s"
probeText += f"{shock_speed_loc} - {gData["shock-speed"][shock_speed_loc][0]/1000:.2f} $\\pm${gData["shock-speed"][shock_speed_loc][1]/1000:.2f} [{gData["shock-speed"][shock_speed_loc][1]/gData["shock-speed"][shock_speed_loc][0]*100:.2f}%] km/s"
graphData["plots"].append({
"type": "text",
"text": f"Measured Shock Speeds{probeText}",
"align": ("top", "right"),
"x": 0.94, "y": 0.94
"x": 0.94 if len(gData["info"]["probe-info"]["locations"]) < 3 else 0.90,
"y": 0.94
})
if doLimits and len(lims) > 1:
@ -420,7 +471,7 @@ print("Graphing Data")
for shot in loaded_data:
#if shot != loaded_data[-2]: continue
genGraph(data[shot], showPlot=False)
genGraph(data[shot], showPlot=False, doLimits=False, forcePlots=True)
#genGraph(data[shot], showPlot=False, doLimits=False, forcePlots=True)
# This forces matplotlib to hang until I tell it to close all windows
pltKeyClose()

View File

@ -35,3 +35,37 @@ volt-scale:
trigbox: 0.001 #V / mV
trigbox_delay: 0.001 #V / mV
uncertainties:
probe-dist:
sd1: 0.002 # +/- m
sd2: 0.002 # +/- m
sd3: 0.002 # +/- m
st1: 0.002 # +/- m
st2: 0.002 # +/- m
st3: 0.002 # +/- m
at1: 0.002 # +/- m
at2: 0.002 # +/- m
at3: 0.002 # +/- m
at4: 0.002 # +/- m
at5: 0.002 # +/- m
at6: 0.002 # +/- m
# [TODO] This could be better
st1-g1: 0.002 # +/- m
st1-g2: 0.002 # +/- m
st2-g1: 0.002 # +/- m
st2-g2: 0.002 # +/- m
st3-g1: 0.002 # +/- m
st3-g2: 0.002 # +/- m
time:
x2-daq: 0