添加相对饱和度计算

.vscode/settings.json

@@ -69,6 +69,7 @@
         "bak": "c",
         "any": "cpp",
         "unordered_set": "cpp",
-        "csignal": "cpp"
+        "csignal": "cpp",
+        "zh2o2_computer.h": "c"
     }
 }

CMakeLists.txt

@@ -13,6 +13,8 @@ set_property(GLOBAL PROPERTY RULE_LAUNCH_COMPILE ccache)
 set_property(GLOBAL PROPERTY RULE_LAUNCH_LINK ccache)
 
 # 设置通用编译选项
+
+set(C_CPP_FLAGS "${C_CPP_FLAGS} -rdynamic")
 set(C_CPP_FLAGS "${C_CPP_FLAGS} -Wno-format-overflow")
 set(C_CPP_FLAGS "${C_CPP_FLAGS} -Wno-unused-local-typedefs")
 set(C_CPP_FLAGS "${C_CPP_FLAGS} -Wno-unused-but-set-variable")

app_protocols/appexception/appexception.hpp

@@ -53,7 +53,10 @@ class appexception : public std::exception {
     return ret;
   }
 
-  appexception() {}
+  appexception() {
+    this->traceinfo = printStackTrace();
+    whatstr         =     this->traceinfo;
+  }
 
   appexception(ThisClass &thisclass, string func, int32_t ecode, const string &description) {
     this->description = description;

appsrc/appbase/utils/h2o2_computer/zh2o2_computer.c

@@ -0,0 +1,207 @@
+#include "zh2o2_computer.h"
+
+#include <stdio.h>
+/**
+ * @brief 水蒸气饱和压力(Pws)
+ *
+ * @param pressure
+ * @return float
+ *
+ * Demo:
+ *  waterVaporSaturationPressure(20+273.15) = 2339.25
+ *
+ */
+static double fn05_Pws(double Tk) {
+  // pws(Tk) = a0* exp(a1/Tk + a2 + a3*Tk + a4*Tk^2 + a5*ln(a6*Tk))
+
+  const double a0 = 1.000000;
+  const double a1 = -6.0969385 * 1000;
+  const double a2 = 2.12409642 * 10;
+  const double a3 = -2.711193 * 0.01;
+  const double a4 = 1.673952 * 0.00001;
+  const double a5 = 2.433502;
+  const double a6 = 1.000000;
+
+  double pws = a0 * exp(a1 / Tk + a2 + a3 * Tk + a4 * Tk * Tk + a5 * log(a6 * Tk));
+
+  return pws;
+}
+/**
+ * @brief
+ */
+
+/**
+ * @brief H2O2蒸气饱和压力Phs
+ *
+ * @param Tk
+ * @return double
+ *
+ * Formula:
+ *          Phs(Tk) = 10^(a0+a1/Tk+a2*log10(Tk)+a3*Tk) / (7.5006*0.001)
+ *
+ * Demo   :
+ *          h2o2VaporSaturationPressure(20+273.15) = 181.295
+ *
+ */
+static double fn06_Phs(double Tk) {
+  const double a0 = 44.576;
+  const double a1 = -4025;
+  const double a2 = -12.996;
+  const double a3 = 0.004605;
+
+  double phs = pow(10, a0 + a1 / Tk + a2 * log10(Tk) + a3 * Tk) / (7.5006 * 0.001);
+  return phs;
+}
+
+/**
+ * @brief fn01_Ph
+ *
+ * @param p
+ * @param pws
+ * @return double
+ *
+ */
+static double fn01_Ph(double ppm, double atmospheric_pressure) { return ppm * atmospheric_pressure / 1000000.0; }
+
+/**
+ * @brief fn10_Yw Eq10(Page7)
+ *
+ * @param Xw
+ * @return double
+ */
+static double fn10_Yw(double Xw, double Tk) {
+  // Yw = exp{ (1-Xw)^2 * (1/RTk) * [B0+B1(1-4Xw)+B2(1-2Xw)*(1-6Xw)] }
+  // R  = Universal gas constant 1.986 cal∙K-1
+
+  double       b0 = -1017 + 0.97 * Tk;
+  const double b1 = 85;
+  const double b2 = 13;
+  const double R  = 1.986;
+
+  double Yw = exp((1 - Xw) * (1 - Xw) * (1 / Tk) * (b0 + b1 * (1 - 4 * Xw) + b2 * (1 - 2 * Xw) * (1 - 6 * Xw)) / R);
+
+  return Yw;
+}
+
+/**
+ * @brief fn11_Yh Eq11(Page7)
+ *
+ * @param Xw
+ * @param Tk
+ * @return double
+ */
+static double fn11_Yh(double Xw, double Tk) {
+  // Yh = exp{ Xw^2/RTk * [B0+B1(3-4Xw)+B2(1-2Xw)*(5-6Xw)] }
+
+  double       b0 = -1017 + 0.97 * Tk;
+  const double b1 = 85;
+  const double b2 = 13;
+  const double R  = 1.986;
+
+  double Yh = exp((Xw * Xw / (Tk * R)) * (b0 + b1 * (3 - 4 * Xw) + b2 * (1 - 2 * Xw) * (5 - 6 * Xw)));
+  return Yh;
+}
+
+static double fn8_Pw(double Pws, double Phs, double Ph, double Yw, double Yh) { return Yw * Pws * (1 - Ph / (Yh * Phs)); }
+static double fn17_Pms(double Pws, double Xw, double Yw) { return Yw * Xw * Pws; }
+
+static void _iteration_pms_and_pw(double xh_avg, double h2o2ppm, double Tk, double AirPressure, double* Pms, double* Pw) {  //
+  double Xw = 1 - xh_avg;
+  double Yw = fn10_Yw(Xw, Tk);
+  double Yh = fn11_Yh(Xw, Tk);
+
+  double Pws = fn05_Pws(Tk);
+  double Phs = fn06_Phs(Tk);
+
+  double Ph = fn01_Ph(h2o2ppm, AirPressure);
+
+  *Pw  = fn8_Pw(Pws, Phs, Ph, Yw, Yh);
+  *Pms = fn17_Pms(Pws, Xw, Yw);
+}
+
+static double t2k(double t) { return t + 273.15; }
+
+/*******************************************************************************
+ *                                   Extern                                    *
+ *******************************************************************************/
+
+double zh2o2_compute_pms(double ppm, double Tk, double AirPressure) {
+  double Xmin = 0;
+  double Xmax = 1;
+  double Pms, Pw;
+
+  int iteration_times = 0;
+  while (true) {
+    iteration_times++;
+
+    if (iteration_times > 100) {
+      printf("%s:%d warninging......,iteration_times > 100\n", __FILE__, __LINE__);
+      break;
+    }
+
+    double xh_avg = (Xmin + Xmax) / 2;
+    if (fabs(Xmax - Xmin) < 0.00001) {
+      break;
+    }
+
+    _iteration_pms_and_pw(xh_avg, ppm, Tk, AirPressure, &Pms, &Pw);
+    // printf("%d xh_avg=%lf,Pms=%lf,Pw=%lf\n", iteration_times, xh_avg, Pms, Pw);
+
+    if (Pw > Pms) {
+      Xmax = xh_avg;
+    } else {
+      Xmin = xh_avg;
+    }
+  }
+
+  return Pms;
+}
+
+double zh2o2_compute_rs(double ppm, double Tk, double rh100, double AirPressure) {
+  double Pms   = zh2o2_compute_pms(ppm, Tk, AirPressure);
+  double Pw    = fn05_Pws(Tk) * rh100 / 100;
+  float  rs100 = Pw / Pms * 100.0;
+  if (rs100 > 100) rs100 = 100;
+  if (rs100 < 0) rs100 = 0;
+  return rs100;
+}
+
+double zh2o2_t2k(double temperatureC) { return temperatureC + 273.15; }
+
+#define ZCHECK_EQ(val0, val1, permision)                                                                                    \
+  {                                                                                                                         \
+    double val0_ = val0;                                                                                                    \
+    double val1_ = val1;                                                                                                    \
+    if (fabs(val0_ - val1_) > permision) {                                                                                  \
+      printf("%s:%d ZCHECK_EQ(%s,%s,%s) failed: %lf != %lf\n", __FILE__, __LINE__, #val0, #val1, #permision, val0_, val1_); \
+      return -1;                                                                                                            \
+    }                                                                                                                       \
+  }
+
+int zh2o2_computer_self_test() {
+  ZCHECK_EQ(fn05_Pws(t2k(20)), 2339.25, 0.01);
+  ZCHECK_EQ(fn06_Phs(t2k(20)), 181.295, 0.01);
+  ZCHECK_EQ(fn01_Ph(500, 101325), 50.6625, 0.01);
+  ZCHECK_EQ(fn10_Yw(0.5, t2k(20)), 0.703911, 0.01);
+  ZCHECK_EQ(fn11_Yh(0.5, t2k(20)), 0.757218, 0.01);
+
+  double Pms, Pw;
+  _iteration_pms_and_pw(0.5, 500, t2k(20), 101325, &Pms, &Pw);
+  ZCHECK_EQ(Pw, 1038.94585, 0.01);
+  ZCHECK_EQ(Pms, 823.311334, 0.01);
+
+  ZCHECK_EQ(zh2o2_compute_pms(500, t2k(20), 101325), 1063.920083, 0.01);
+  ZCHECK_EQ(zh2o2_compute_rs(500, t2k(20), 42, 101325), 92.345718, 0.01);
+
+  float temperatureC = 26;
+  float AirPressure  = 101325;
+  float rh           = 50;
+  printf("T  AirPressure RH  H2O2(ppm)  RS\n");
+  for (size_t i = 0; i < 500; i += 1) {
+    float rs = zh2o2_compute_rs(i, zh2o2_t2k(temperatureC), rh, AirPressure);
+    printf("%d %d %d %d %d\n", (int) temperatureC, (int)AirPressure,  (int)rh, (int) i,  (int)rs);
+  }
+
+  printf("Test passed\n");
+  return 0;
+}

appsrc/appbase/utils/h2o2_computer/zh2o2_computer.h

@@ -0,0 +1,46 @@
+#pragma once
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <math.h>
+#include <stdbool.h>
+#include <stdint.h>
+
+//https://iflytop1.feishu.cn/wiki/QvXow7g6aiI1kskUVwQcRL1BnXd
+
+#define AIR_PRESSURE 101325.0
+
+/**
+ * @brief 计算 water vapor saturation pressure  in an H2O2 vapor environment
+ *
+ * @param ppm 过氧化氢浓度
+ * @param Tk  温度(K)
+ * @param AirPressure 大气压
+ * @return double
+ */
+double zh2o2_compute_pms(double ppm, double Tk, double AirPressure);
+/**
+ * @brief 计算 过氧气化氢相对饱和度(RS)
+ *
+ * @param ppm 过氧化氢浓度
+ * @param Tk  温度(K)
+ * @param rh100 相对湿度(RH%)
+ * @param AirPressure 大气压
+ * @return double
+ */
+double zh2o2_compute_rs(double ppm, double Tk, double rh100, double AirPressure);
+/**
+ * @brief 计算 摄氏度转开尔文
+ *
+ * @param temperatureC
+ * @return double
+ */
+double zh2o2_t2k(double temperatureC);
+
+int zh2o2_computer_self_test();
+
+#ifdef __cplusplus
+}
+#endif

appsrc/appcomponents/canchannel/protocol_utils.hpp

@@ -0,0 +1,5 @@
+#pragma once
+#include "app_protocols/transmit_disfection_protocol/transmit_disfection_protocol.hpp"
+
+namespace iflytop { 
+}

appsrc/service/hardware/base/h2o2_sensor_data_mgr.cpp

@@ -1,4 +1,6 @@
 #include "h2o2_sensor_data_mgr.hpp"
+
+#include "appbase/utils/h2o2_computer/zh2o2_computer.h"
 using namespace iflytop;
 
 bool H2O2SensorDataMgr::isDataExpired(uint8_t sensorId) { return isDataExpired(getSensorDataCache(sensorId)); }
@@ -83,7 +85,10 @@ void H2O2SensorDataMgr::updateH2o2SensorData(uint8_t sensorId, report_h2o2_data_
     float h2o2 = h2o2data->h2o2;
     float temp = h2o2data->temp / 10.0;
     float rh   = h2o2data->rh / 10.0;
-    float rs   = h2o2data->rs / 10.0;
+    float rs   = 0;
+    if (h2o2data->rs == 0) {
+      rs = zh2o2_compute_rs(h2o2, zh2o2_t2k(temp), rh, AIR_PRESSURE);
+    }
     if (h2o2 < 30) {
       h2o2 = 0;
     }
@@ -123,7 +128,7 @@ void H2O2SensorDataMgr::statisticsSensorData() {
       continue;
     }
 
-    if(m_minH2O2 == -1) m_minH2O2 = sensor_data->h2o2;
+    if (m_minH2O2 == -1) m_minH2O2 = sensor_data->h2o2;
 
     if (sensor_data->h2o2 < m_minH2O2) {
       m_minH2O2 = sensor_data->h2o2;

test/zh2o2_computer_selftest.sh

@@ -0,0 +1,7 @@
+#!/bin/bash
+set -e
+echo "extern int zh2o2_computer_self_test();" > /tmp/tmpmain.c
+echo "int main(int argc, char* const argv[]) { zh2o2_computer_self_test(); }" >> /tmp/tmpmain.c
+gcc ../appsrc/appbase/utils/h2o2_computer/zh2o2_computer.c /tmp/tmpmain.c -lm -o selftest.out
+./selftest.out
+rm selftest.out