"""气体流量开环扫点实验(不参与闭环反馈,也不被 main.py 导入)。 本模块可独立运行,用于阀门开度—流量特性的开环标定/扫点。程序会: 1. 生成 0%~100%(默认 5% 间隔)共 21 个开度值并随机打乱; 2. 逐个开度:写入电机行程 → 以固定周期采集控流阀前后 4 个传感器 (前压力、前流量、后压力、后流量)→ 等待“控流阀后流量”进入稳态后 切换到下一个开度; 3. 全程逐周期写 CSV,结束后在单个窗口内输出三张子图(前后压力、前后流量、 绝对压力比),保存到工程目录下的 open_loop_data/。 示例: python open_loop.py python open_loop.py --seed 42 # 复现同一种打乱顺序 python open_loop.py --min 0 --max 100 --step 5 量程、行程端值、稳态窗口等沿用 config.py;本文件顶部的常量可在不修改 config.py 的前提下覆盖实验相关参数(尤其是四个传感器通道映射)。 """ import argparse import csv import logging import math import random import statistics import time from collections import deque from pathlib import Path import config from flow_control import opening_to_motor_position # --------------------------------------------------------------------------- # 实验参数(可在本文件修改,不影响 config.py) # --------------------------------------------------------------------------- # 四个传感器的 AI 通道(0 起始,MT2-AM8 AI0~AI3)。 # “后”= 控流阀出口侧(靠近流量计/负载一侧),沿用 config 的原始通道; # “前”= 控流阀入口侧(靠近减压阀一侧),为本次新增传感器,请按实际接线确认。 PRESSURE_AFTER_CHANNEL = config.PRESSURE_INPUT_CHANNEL # 后压力,默认 0 FLOW_AFTER_CHANNEL = config.FLOW_INPUT_CHANNEL # 后流量,默认 1 PRESSURE_BEFORE_CHANNEL = None # 前压力(新增,需确认) FLOW_BEFORE_CHANNEL = None # 前流量(新增,需确认) MOTOR_CHANNEL = config.MOTOR_OUTPUT_CHANNEL # AO 通道 # 开度扫描:0%~100%,默认 5% 间隔,共 21 个开度值,随机打乱。 OPENING_MIN_PCT = 0.0 OPENING_MAX_PCT = 100.0 OPENING_STEP_PCT = 5.0 RANDOM_SEED = 42 # 设为整数可复现同一打乱顺序 # 采样与稳态判定 SAMPLE_PERIOD_S = 0.1 # 采样周期(秒) STEADY_STATE_WINDOW_S = config.STEADY_STATE_WINDOW_S # 滑窗长度,复用 config STEADY_STATE_STD_PCT = config.STEADY_STATE_STD_PCT # 相对标准差阈值(%) MAX_WAIT_S = 60.0 # 单个开度的最长等待时间 MIN_ABS_STD_SLM = 0.05 # 近零流量时的绝对标准差下限 MAX_CONSECUTIVE_READ_FAILURES = 5 # 后流量连续读取失败判定设备断开 # 安全阈值(沿用 config 的闭环安全值;0% 开度时近零噪声若略为负,可放宽流量下限) SAFETY_FLOW_MIN_SLM = config.FLOW_VALID_MIN_SLM SAFETY_FLOW_MAX_SLM = config.FLOW_VALID_MAX_SLM SAFETY_MAX_PRESSURE_KPA = config.MAX_PRESSURE_KPA # 输出 OUTPUT_DIRECTORY = "open_loop_data" # 相对工程目录 PLOT_DPI = config.PLOT_DPI ANNOTATE_STEPS = True # 曲线上标注开度步 # 绝对压力 = 表压读数 + 大气压;临界压力比用于判断是否出现壅塞(choked)流动。 ATMOSPHERIC_PRESSURE_KPA = 101.325 # 标准大气压(kPa) CRITICAL_PRESSURE_RATIO = 0.528 # 临界压力比(阀后绝对压力 / 阀前绝对压力) CSV_FIELDS = ( "time_s", "step_index", "opening_pct", "motor_position", "flow_before_slm", "flow_after_slm", "pressure_before_kpa", "pressure_after_kpa", "P_abs_ratio", "is_ratio_smaller_than_0.528", ) # --------------------------------------------------------------------------- # 稳态判定 # --------------------------------------------------------------------------- class SteadyStateDetector: """滑动窗口稳态判定,复用 PerformanceReporter 的滑窗标准差思路。 闭环有目标流量,所以用“均值落在目标 ±band”+“标准差低于阈值”两个条件; 开环没有目标,因此只保留“滑窗相对标准差足够小”这一条,表示信号不再 明显变化: std <= max(STD_PCT% * |mean|, MIN_ABS_STD_SLM) 近零流量时以绝对下限兜底,避免除以零。 """ def __init__(self, window_s, std_pct, sample_period_s=0.1, min_abs_std_slm=0.05): self.window_s = float(window_s) self.std_pct = float(std_pct) self.sample_period_s = float(sample_period_s) self.min_abs_std_slm = float(min_abs_std_slm) self._window = deque() # 元素为 (elapsed_s, flow_slm) def reset(self): self._window.clear() def observe(self, elapsed_s, flow_slm): """喂入一个样本;达到稳态返回 True,否则返回 False。""" if flow_slm is None or not math.isfinite(float(flow_slm)): return False self._window.append((float(elapsed_s), float(flow_slm))) cutoff = float(elapsed_s) - self.window_s while self._window and self._window[0][0] < cutoff: self._window.popleft() if len(self._window) < 2: return False # 滑窗覆盖时长需至少达到 window_s。真实时间戳带抖动,若严格比较 # span < window_s,几乎总差一个采样间隙而永远无法判定;放宽一个 # 采样周期的余量。 min_span = max(0.0, self.window_s - self.sample_period_s) if self._window[-1][0] - self._window[0][0] < min_span: return False flows = [item[1] for item in self._window] mean = statistics.mean(flows) std = statistics.pstdev(flows) threshold = max(self.std_pct / 100.0 * abs(mean), self.min_abs_std_slm) return std <= threshold # --------------------------------------------------------------------------- # 数据记录与绘图 # --------------------------------------------------------------------------- class OpenLoopRecorder: """逐周期写 CSV,结束时输出三张子图(前后压力、前后流量、绝对压力比)。""" def __init__(self, output_directory, *, plot_dpi=160): output_dir = Path(output_directory) output_dir.mkdir(parents=True, exist_ok=True) run_id = time.strftime("%Y%m%d_%H%M%S") self.csv_path = output_dir / f"open_loop_{run_id}.csv" self.image_path = output_dir / f"open_loop_{run_id}.png" self.plot_dpi = int(plot_dpi) self._file = self.csv_path.open("w", newline="", encoding="utf-8-sig") self._writer = csv.DictWriter(self._file, fieldnames=CSV_FIELDS) self._writer.writeheader() self._file.flush() self.sample_count = 0 self.steps = [] # 每个开度步的 (time_s, opening_pct) self._last_step_index = None self._closed = False def record( self, *, time_s, step_index, opening_pct, motor_position, flow_before_slm, flow_after_slm, pressure_before_kpa, pressure_after_kpa, ): if self._closed: raise RuntimeError("记录器已关闭") p_abs_ratio = self._absolute_pressure_ratio( pressure_before_kpa, pressure_after_kpa ) if p_abs_ratio is None: p_abs_ratio_text = "" is_smaller_text = "" else: p_abs_ratio_text = f"{p_abs_ratio:.6f}" is_smaller_text = ( "Y" if p_abs_ratio < CRITICAL_PRESSURE_RATIO else "N" ) self._writer.writerow( { "time_s": self._fmt(time_s), "step_index": int(step_index), "opening_pct": self._fmt(opening_pct), "motor_position": self._fmt(motor_position), "flow_before_slm": self._fmt(flow_before_slm), "flow_after_slm": self._fmt(flow_after_slm), "pressure_before_kpa": self._fmt(pressure_before_kpa), "pressure_after_kpa": self._fmt(pressure_after_kpa), "P_abs_ratio": p_abs_ratio_text, "is_ratio_smaller_than_0.528": is_smaller_text, } ) self._file.flush() self.sample_count += 1 if step_index != self._last_step_index: self.steps.append((float(time_s), float(opening_pct))) self._last_step_index = step_index @staticmethod def _absolute_pressure_ratio(pressure_before_kpa, pressure_after_kpa): """阀后绝对压力 / 阀前绝对压力;读数缺失或分母非法时返回 None。""" if pressure_before_kpa is None or pressure_after_kpa is None: return None try: before_abs = float(pressure_before_kpa) + ATMOSPHERIC_PRESSURE_KPA after_abs = float(pressure_after_kpa) + ATMOSPHERIC_PRESSURE_KPA except (TypeError, ValueError): return None if before_abs <= 0: return None ratio = after_abs / before_abs return ratio if math.isfinite(ratio) else None def finalize(self): """关闭 CSV 并生成 PNG;无采样点时只保留 CSV。""" self.close() if self.sample_count == 0: return None self._create_plot() return self.image_path def close(self): if not self._closed: self._file.flush() self._file.close() self._closed = True def _create_plot(self): # 无界面后端,保证终端/无显示器环境也能保存图片。 import matplotlib matplotlib.use("Agg") import matplotlib.pyplot as plt time_s = [] flow_before = [] flow_after = [] pressure_before = [] pressure_after = [] p_abs_ratio = [] with self.csv_path.open("r", newline="", encoding="utf-8-sig") as file: for row in csv.DictReader(file): time_s.append(float(row["time_s"])) flow_before.append(self._float_or_nan(row["flow_before_slm"])) flow_after.append(self._float_or_nan(row["flow_after_slm"])) pressure_before.append( self._float_or_nan(row["pressure_before_kpa"]) ) pressure_after.append( self._float_or_nan(row["pressure_after_kpa"]) ) p_abs_ratio.append(self._float_or_nan(row["P_abs_ratio"])) figure, axes = plt.subplots( 3, 1, figsize=(12, 12), sharex=True, constrained_layout=True, ) figure.suptitle("Open-loop Valve Sweep", fontsize=15) self._plot_series( axes[0], time_s, pressure_before, "Pressure before valve", "#1565C0", ) self._plot_series( axes[0], time_s, pressure_after, "Pressure after valve", "#D84315", ) axes[0].set_ylabel("Pressure (kPa)") axes[0].set_title("Pressure across control valve") axes[0].legend(loc="best") self._plot_series( axes[1], time_s, flow_before, "Flow before valve", "#2E7D32", ) self._plot_series( axes[1], time_s, flow_after, "Flow after valve", "#6A1B9A", ) axes[1].set_ylabel("Flow (SLM)") axes[1].set_title("Flow across control valve") axes[1].legend(loc="best") self._plot_series( axes[2], time_s, p_abs_ratio, "P_abs_ratio (after/before)", "#EF6C00", ) axes[2].axhline( CRITICAL_PRESSURE_RATIO, color="#B71C1C", linestyle="--", linewidth=1.2, label=f"Critical ratio {CRITICAL_PRESSURE_RATIO:.3f}", ) axes[2].set_ylabel("P_abs_ratio") axes[2].set_xlabel("Time (s)") axes[2].set_title("Absolute pressure ratio across control valve") axes[2].legend(loc="best") if ANNOTATE_STEPS and self.steps: y_top = axes[1].get_ylim()[1] for step_time, opening in self.steps: for axis in axes: axis.axvline( step_time, color="0.5", linestyle="--", linewidth=0.7, alpha=0.6, ) axes[1].text( step_time, y_top, f"{opening:.0f}%", rotation=90, fontsize=7, ha="left", va="bottom", color="0.3", ) for axis in axes: axis.grid(True, alpha=0.3, linestyle="--") axis.margins(x=0) figure.savefig(self.image_path, dpi=self.plot_dpi, bbox_inches="tight") plt.close(figure) @staticmethod def _plot_series(axis, time_s, series, label, color): if any(math.isfinite(value) for value in series): axis.plot(time_s, series, color=color, linewidth=1.4, label=label) else: axis.text( 0.5, 0.5, f"No data: {label}", transform=axis.transAxes, ha="center", va="center", color="gray", ) @staticmethod def _fmt(value): if value is None: return "" number = float(value) return "" if not math.isfinite(number) else f"{number:.6f}" @staticmethod def _float_or_nan(value): if value in (None, ""): return math.nan try: return float(value) except (TypeError, ValueError): return math.nan # --------------------------------------------------------------------------- # 硬件读写辅助 # --------------------------------------------------------------------------- def _safe_read(channel, read_fn): """读取单路传感器;通道为 None 或读取失败时返回 None,不抛异常。""" if channel is None: return None try: return read_fn(channel) except Exception: return None def read_sensors(hardware): """读取控流阀前后四个传感器,返回 (前压, 前流量, 后压, 后流量)。""" pressure_before = _safe_read(PRESSURE_BEFORE_CHANNEL, hardware.get_pressure) flow_before = _safe_read(FLOW_BEFORE_CHANNEL, hardware.get_flow) pressure_after = _safe_read(PRESSURE_AFTER_CHANNEL, hardware.get_pressure) flow_after = _safe_read(FLOW_AFTER_CHANNEL, hardware.get_flow) return pressure_before, flow_before, pressure_after, flow_after def check_safety(pressure_before, flow_before, pressure_after, flow_after): """越限检查,返回故障描述字符串;正常返回 None。 沿用闭环 FlowControlLoop 的安全阈值:流量超出 [SAFETY_FLOW_MIN_SLM, SAFETY_FLOW_MAX_SLM],或任一压力超过 SAFETY_MAX_PRESSURE_KPA 时判定为故障(顺带识别未接线通道的负读数)。 """ if ( flow_after is not None and not SAFETY_FLOW_MIN_SLM <= flow_after <= SAFETY_FLOW_MAX_SLM ): return f"FLOW_AFTER_OUT_OF_RANGE({flow_after:.3f} SLM)" if ( flow_before is not None and not SAFETY_FLOW_MIN_SLM <= flow_before <= SAFETY_FLOW_MAX_SLM ): return f"FLOW_BEFORE_OUT_OF_RANGE({flow_before:.3f} SLM)" if pressure_after is not None and pressure_after > SAFETY_MAX_PRESSURE_KPA: return f"PRESSURE_AFTER_OVER_LIMIT({pressure_after:.3f} kPa)" if pressure_before is not None and pressure_before > SAFETY_MAX_PRESSURE_KPA: return f"PRESSURE_BEFORE_OVER_LIMIT({pressure_before:.3f} kPa)" return None def set_opening(hardware, opening_pct): """把 0~100% 开度换算成反向电机行程并下发给 MT2-AM8。""" position = opening_to_motor_position( opening_pct, config.MOTOR_OPEN_POSITION, config.MOTOR_CLOSED_POSITION, ) return bool( hardware.set_motor_position(position, channel=MOTOR_CHANNEL) ), position def generate_openings(min_pct, max_pct, step_pct, seed): """生成 [min, max] 内 step 间隔的开度序列并随机打乱。""" count = int(round((max_pct - min_pct) / step_pct)) + 1 openings = [round(min_pct + i * step_pct, 6) for i in range(count)] # step 不能整除区间时,自动补上 max 端点,保证上下限都被覆盖。 if not openings or abs(openings[-1] - max_pct) > 1e-6: openings.append(float(max_pct)) rng = random.Random(seed) rng.shuffle(openings) return openings # --------------------------------------------------------------------------- # 日志 # --------------------------------------------------------------------------- def build_logger(): log_dir = Path(__file__).resolve().parent / config.LOG_DIRECTORY log_dir.mkdir(parents=True, exist_ok=True) log_path = log_dir / "open_loop.log" logger = logging.getLogger("open_loop") logger.setLevel(logging.INFO) logger.handlers.clear() formatter = logging.Formatter( "%(asctime)s.%(msecs)03d %(levelname)s %(message)s", datefmt="%Y-%m-%d %H:%M:%S", ) file_handler = logging.FileHandler(log_path, encoding="utf-8") file_handler.setFormatter(formatter) logger.addHandler(file_handler) console_handler = logging.StreamHandler() console_handler.setFormatter(formatter) logger.addHandler(console_handler) return logger # --------------------------------------------------------------------------- # 主流程 # --------------------------------------------------------------------------- def parse_args(argv=None): parser = argparse.ArgumentParser(description="MT2-AM8 气体流量开环扫点实验") parser.add_argument("--min", type=float, default=OPENING_MIN_PCT, help=f"扫描起始开度(%%),默认 {OPENING_MIN_PCT:.0f}") parser.add_argument("--max", type=float, default=OPENING_MAX_PCT, help=f"扫描结束开度(%%),默认 {OPENING_MAX_PCT:.0f}") parser.add_argument("--step", type=float, default=OPENING_STEP_PCT, help=f"开度间隔(%%),默认 {OPENING_STEP_PCT:.0f}") parser.add_argument("--seed", type=int, default=RANDOM_SEED, help="随机种子,默认 None(每次打乱顺序不同)") args = parser.parse_args(argv) if args.step <= 0: parser.error("--step 必须大于 0") if not (0 <= args.min < args.max <= 100): parser.error("必须满足 0 <= --min < --max <= 100") return args def run(args): config.validate_config() logger = build_logger() try: from PcControl import MT2AM8Client except ModuleNotFoundError as exc: if exc.name and exc.name.startswith("pymodbus"): logger.critical( "缺少 PcControl.py 所需的 pymodbus;请在运行本工程的 Python " "环境中安装与现有硬件代码兼容的 pymodbus 版本" ) return 6 raise hardware = MT2AM8Client( host=config.MT2AM8_HOST, port=config.MT2AM8_PORT, slave_id=config.MT2AM8_SLAVE_ID, pressure_range=config.PRESSURE_RANGE_KPA, flow_range=config.FLOW_METER_RANGE_SLM, ) openings = generate_openings( args.min, args.max, args.step, args.seed ) logger.info( "开度扫描顺序(已打乱): %s", ", ".join(f"{o:.1f}%" for o in openings), ) connected = False recorder = None exit_code = 0 try: logger.info("正在连接 MT2-AM8 %s:%s", config.MT2AM8_HOST, config.MT2AM8_PORT) connected = bool(hardware.connect()) if not connected: logger.critical("无法连接 MT2-AM8,实验中止") exit_code = 2 return exit_code # 启动前先关阀。 ok, _ = set_opening(hardware, 0.0) if not ok: logger.critical("启动前关阀失败,实验中止") exit_code = 2 return exit_code recorder = OpenLoopRecorder( Path(__file__).resolve().parent / OUTPUT_DIRECTORY, plot_dpi=PLOT_DPI, ) logger.info("逐周期数据将保存到 %s", recorder.csv_path) detector = SteadyStateDetector( STEADY_STATE_WINDOW_S, STEADY_STATE_STD_PCT, sample_period_s=SAMPLE_PERIOD_S, min_abs_std_slm=MIN_ABS_STD_SLM, ) run_start = time.perf_counter() read_failures = 0 for step_index, opening in enumerate(openings): ok, position = set_opening(hardware, opening) if not ok: logger.critical( "写入开度 %.1f%% 对应行程 %.1f 失败,实验中止", opening, position ) exit_code = 2 break logger.info( "步 %d/%d:设置开度 %.1f%%(行程 %.1f),等待后流量稳定…", step_index + 1, len(openings), opening, position, ) detector.reset() step_start = time.perf_counter() next_sample = step_start stable = False while True: now = time.perf_counter() elapsed_s = now - run_start pressure_before, flow_before, pressure_after, flow_after = ( read_sensors(hardware) ) recorder.record( time_s=elapsed_s, step_index=step_index, opening_pct=opening, motor_position=position, flow_before_slm=flow_before, flow_after_slm=flow_after, pressure_before_kpa=pressure_before, pressure_after_kpa=pressure_after, ) fault = check_safety( pressure_before, flow_before, pressure_after, flow_after ) if fault is not None: logger.critical("安全阈值触发:%s,实验中止并关阀", fault) exit_code = 2 break if flow_after is None: read_failures += 1 else: read_failures = 0 if read_failures >= MAX_CONSECUTIVE_READ_FAILURES: logger.critical( "控流阀后流量连续 %d 次读取失败,判定设备断开,实验中止", read_failures, ) exit_code = 2 break if detector.observe(elapsed_s, flow_after): stable = True settle_s = now - step_start flow_text = "--" if flow_after is None else f"{flow_after:.3f}" logger.info( "步 %d 完成:开度 %.1f%% 后流量稳定 %s SLM,耗时 %.2f s", step_index + 1, opening, flow_text, settle_s, ) break if now - step_start >= MAX_WAIT_S: logger.warning( "开度 %.1f%% 在 %.1f s 内未稳定,强制进入下一开度", opening, MAX_WAIT_S, ) break next_sample += SAMPLE_PERIOD_S sleep_s = next_sample - time.perf_counter() if sleep_s > 0: time.sleep(sleep_s) else: next_sample = time.perf_counter() if exit_code != 0: break except KeyboardInterrupt: logger.info("收到 Ctrl+C,正在停止实验") except Exception: exit_code = 3 logger.exception("未处理异常,实验进入收尾") finally: if connected: try: close_ok, _ = set_opening(hardware, 0.0) if not close_ok: exit_code = max(exit_code, 4) logger.critical("收尾安全关阀失败,请立即人工检查") except Exception: exit_code = max(exit_code, 4) logger.exception("收尾安全关阀时发生异常") try: hardware.disconnect() except Exception: exit_code = max(exit_code, 5) logger.exception("断开 MT2-AM8 时发生异常") if recorder is not None: try: image_path = recorder.finalize() logger.info("开环数据已保存:%s", recorder.csv_path) if image_path is not None: logger.info("曲线图已保存:%s", image_path) else: logger.warning("本次运行没有采样点,因此未生成曲线图") except Exception: recorder.close() exit_code = max(exit_code, 7) logger.exception("保存曲线图失败;CSV 数据仍保留") logger.info("程序结束,退出码=%d", exit_code) return exit_code def main(argv=None): args = parse_args(argv) return run(args) if __name__ == "__main__": import sys sys.exit(main())