// ignore_for_file: constant_identifier_names import 'dart:math' as math; import 'package:fis_measure/configs/patient.dart'; import 'package:fis_measure/interfaces/date_types/point.dart'; import 'package:fis_measure/interfaces/enums/species.dart'; class GeneralFormulas { GeneralFormulas._(); static const double VolumeCofficient = math.pi / 6.0; static final IGeneralFormulaStrategy _singleton = GlobalPatientConfig.speciesType == SpeciesType.mouse ? BaseGeneralFormulas() : AnimalsGeneralFormulas(); static double volumeTwoLine(double d1, double d2, [bool useBigValue = true]) => _singleton.volumeTwoLine(d1, d2, useBigValue); static double volumeWithCoefficient( double d1, double d2, double d3, double coefficient) => _singleton.volumeWithCoefficient(d1, d2, d3, coefficient); static double countStenosis(double d1, double d2) => _singleton.countStenosis(d1, d2); static bool doubleAlmostEquals(double num1, double num2, [double precision = 0.000001]) => _singleton.doubleAlmostEquals(num1, num2, precision); static double distance2Line(DPoint pt1, DPoint pt2, DPoint pt) => _singleton.distance2Line(pt1, pt2, pt); static double countSlope(DPoint p1, DPoint p2) => _singleton.countSlope(p1, p2); static double countEnvelopeTime(List points) => _singleton.countEnvelopeTime(points); static double countPressure(double v) => _singleton.countPressure(v); static double countRI(double ps, double mdEd) => _singleton.countRI(ps, mdEd); static double countMaxPG(double v1, double v2) => _singleton.countMaxPG(v1, v2); static double countHR(double timeSpan) => _singleton.countHR(timeSpan); static List countVTI(List tracePoints) => _singleton.countVTI(tracePoints); static double svDiam(double diam, double vti) => _singleton.svDiam(diam, vti); static double csa(double diam) => _singleton.csa(diam); static double dpDtRatio(DPoint p1, DPoint p2) => _singleton.dpDtRatio(p1, p2); static double maxPG(double v1, double v2) => _singleton.maxPG(v1, v2); static double area(double d1, double d2) => _singleton.area(d1, d2); static double perimeter(double d1, double d2) => _singleton.perimeter(d1, d2); static double flowVolTAMAX( double flowArea, double tamean, double coefficient) => _singleton.flowVolTAMAX(flowArea, tamean, coefficient); static double pi(double ps, double mdEd, double tamax) => _singleton.pi(ps, mdEd, tamax); static double medianVelocity(double ps, double ed) => _singleton.medianVelocity(ps, ed); static double volume(double d) => _singleton.volume(d); } abstract class IGeneralFormulaStrategy { double volumeTwoLine(double d1, double d2, [bool useBigValue = true]); double volumeWithCoefficient( double d1, double d2, double d3, double coefficient, ); double countStenosis(double d1, double d2); bool doubleAlmostEquals(double num1, double num2, [double precision = 0.000001]); double distance2Line(DPoint pt1, DPoint pt2, DPoint pt); double countSlope(DPoint p1, DPoint p2); double countEnvelopeTime(List points); double countPressure(double v); double countRI(double ps, double mdEd); double countMaxPG(double v1, double v2); double countHR(double timeSpan); List countVTI( List tracePoints, ); double svDiam(double diam, double vti); double csa(double diam); double dpDtRatio(DPoint p1, DPoint p2); double maxPG(double v1, double v2); double area(double d1, double d2); double flowVolTAMAX(double flowArea, double tamean, double coefficient); double perimeter(double d1, double d2); double volume(double d); /// /// Pulastility Index /// PI (ED):(PS - ED)/TAMAX /// PI (MD):(PS - MD)/TAMAX /// /// Unit cm/s /// Unit cm/s /// Unit cm/s /// Unit none /// 在TCCD下用Vmean替换tamax计算 double pi(double ps, double mdEd, double tamax); /// Unit cm/s /// Unit cm/s /// Unit cm/s double medianVelocity(double ps, double ed); } class BaseGeneralFormulas implements IGeneralFormulaStrategy { /// Volume:1/6 x π static const double VolumeCofficient = math.pi / 6.0; /// 1/6 x π x D1^2 x D2 (D1 > D2) /// /// 1/6 x π x D1 x D2^2 (D1 ≤ D2) /// /// return unit: cm³ @override double volumeTwoLine(double d1, double d2, [bool useBigValue = true]) { double volume = 0; if (!doubleAlmostEquals(d1, 0) && !doubleAlmostEquals(d2, 0)) { if (useBigValue) { if (d1 > d2) { volume = d1 * d1 * d2 * math.pi / 6.0; } else { volume = d1 * d2 * d2 * math.pi / 6.0; } } else { if (d1 > d2) { volume = d2 * d2 * d1 * math.pi / 6.0; } else { volume = d1 * d1 * d2 * math.pi / 6.0; } } } return volume; } /// 带比率系数算体积 /// /// [coefficient] 比率系数 @override double volumeWithCoefficient( double d1, double d2, double d3, double coefficient, ) { double volume = 0; if (!doubleAlmostEquals(d1, 0) && !doubleAlmostEquals(d2, 0) && !doubleAlmostEquals(d3, 0)) { volume = d1 * d2 * d3 * coefficient; } return volume; } /// 狭窄率计算 @override double countStenosis(double d1, double d2) { double residual; double lumen; if (d1 > d2) { residual = d2; lumen = d1; } else { residual = d1; lumen = d2; } return (1.0 - residual / lumen) * 100; } @override bool doubleAlmostEquals(double num1, double num2, [double precision = 0.000001]) { if (num1.isNaN && num2.isNaN) return true; return (num1 - num2).abs() <= precision; } /// 点到线的距离 /// /// [pt1] 线端点1 /// /// [pt2] 线端点2 /// /// [pt] 点 @override double distance2Line(DPoint pt1, DPoint pt2, DPoint pt) { double dis = 0; if (pt1.x == pt2.x) { dis = (pt.x - pt1.x).abs(); return dis; } var lineK = (pt2.y - pt1.y) / (pt2.x - pt1.x); var lineC = (pt2.x * pt1.y - pt1.x * pt2.y) / (pt2.x - pt1.x); dis = (lineK * pt.x - pt.y + lineC).abs() / (math.sqrt(lineK * lineK + 1)); return dis; } ///计算斜率 @override double countSlope(DPoint p1, DPoint p2) { if (doubleAlmostEquals(p2.x, p1.x)) return 0; double vertical = p2.y - p1.y; double time = (p2.x - p1.x).abs(); double slope = vertical / time; return slope; } ///计算包络时间 ///Origin: Vinno.Modules.MeasureModule\Formulas\GeneralFormulas.cs [632:67] @override double countEnvelopeTime(List points) { if (points.length < 2) return 0; double startTime = double.maxFinite; double endTime = -double.maxFinite; for (var point in points) { if (point.x < startTime) { startTime = point.x; } if (point.x > endTime) { endTime = point.x; } } return (endTime - startTime).abs(); } /// /// 4.0 x V^2 /// /// Unit cm/s /// Unit mmHg @override double countPressure(double v) { // The velocity is in cm/s, but it should be m/s in this formula, so divide it by 100. double pg = 4.0 * math.pow(v * 0.01, 2); return pg; } /// /// Resistivity Index 阻力指数 /// RI (ED)=(PS - ED)/PS /// RI (MD)=(PS - MD)/PS /// /// Unit cm/s /// Unit cm/s @override double countRI(double ps, double mdEd) { ps = ps.abs(); mdEd = mdEd.abs(); double ri = (ps - mdEd) / ps; if (ri < 0) { return double.nan; } return ri; } /// /// 4.0 x |V1^2 - V2^2| /// /// Unit cm/s /// Unit cm/s /// Unit mmHg @override double countMaxPG(double v1, double v2) { // The velocity is in cm/s, but it should be m/s in this formula, so divide it by 100. double meanPG = 4.0 * (math.pow(v1 * 0.01, 2) - math.pow(v2 * 0.01, 2)).abs(); return meanPG; } // 计算心率 // TODO: 心脏周期(可配置) DefaultHeartCycle = 2; // Origin: Vinno.Modules.MeasureModule\Primitives\DopplerTraceBase.cs @override double countHR(double timeSpan) { const defaultHeartCycle = 1; double hr = (defaultHeartCycle * 60 / timeSpan).roundToDouble(); return hr; } // return [VTI,tiEnv,timeAveragedVelocity,pv,meanPG,hr]; @override List countVTI( List tracePoints, ) { var tiEnv = 0.0; var timeAveragedVelocity = 0.0; var pv = 0.0; var meanPG = 0.0; if (tracePoints.length < 2) return [0, 0, 0, 0, 0, 0]; int n = tracePoints.length - 1; double dis = 0; double startTime = 0; double endTime = 0; double pgSum = 0; int positiveCount = 0; int negativeCount = 0; for (int i = 0; i < n; i++) { var p1 = tracePoints[i]; var p2 = tracePoints[i + 1]; if (pv.abs() < (p1.y).abs()) { pv = p1.y; } double meanVelocity = (p1.y + p2.y).abs() / 2; double timeSpan = (p1.x - p2.x).abs(); dis += meanVelocity * timeSpan; pgSum += countPressure(meanVelocity) * timeSpan; if (i == 0) { startTime = p1.x; endTime = p2.x; } else { startTime = math.min(startTime, p2.x); endTime = math.max(endTime, p2.x); } if (p1.y < 0) { negativeCount++; } else { positiveCount++; } } if (tracePoints[n].y < 0) { negativeCount++; } else { positiveCount++; } if (positiveCount < negativeCount) { dis = -dis; } if (pv.abs() < (tracePoints[n].y).abs()) { pv = tracePoints[n].y; } tiEnv = (endTime - startTime).abs(); meanPG = (pgSum / tiEnv); timeAveragedVelocity = (dis / tiEnv); var hr = countHR(tiEnv); return [dis, tiEnv, timeAveragedVelocity, pv, meanPG, hr]; } /// SV Diam /// /// SV = 1/4 x π x (Diam)^2 x VTI /// /// [diam] `cm` /// /// [vti] `cm` /// /// return `cm3` @override double svDiam(double diam, double vti) { double sv = math.pi * math.pow(diam, 2) * vti / 4; return sv; } /// CSA /// /// CSA = 1/4 x π x Diam^2 /// /// [diam] `cm` /// /// return `cm²` @override double csa(double diam) { double csa = math.pi * math.pow(diam, 2) / 4; return csa; } /// CW/PW Mode: dp/dt = (4* V2 * V2 - 4 * V1 * V1)/ (T2 - T1) /// /// return `mmHg/s` @override double dpDtRatio(DPoint p1, DPoint p2) { double dp = maxPG(p2.y, p1.y); double dt = (p2.x - p1.x).abs(); double result = dp / dt; return result; } /// 4.0 x |V1^2 - V2^2| /// /// [v1] `cm/s` /// /// [v2] `cm/s` /// /// reutrn `mmHg` @override double maxPG(double v1, double v2) { // The velocity is in cm/s, but it should be m/s in this formula, so divide it by 100. double meanPG = 4.0 * (math.pow(v1 * 0.01, 2) - math.pow(v2 * 0.01, 2)).abs(); return meanPG; } /// 0.25 x π x D1 x D2 /// /// [d1] `cm` /// /// [d2] `cm` /// /// return `cm²` @override double area(double d1, double d2) { double area = 0; if (!doubleAlmostEquals(d1, 0) && !doubleAlmostEquals(d2, 0)) { if ((d1 < 0.001 * d2) || (d2 < 0.001 * d1)) { // Stright line area = 0.0; } else { // Cirle area = 0.25 * d1 * d2 * math.pi; } } return area; } /// 0.25 x π x Diam^2 x TAMAX x Compensation Coefficient, where 0.5 ≤ Compensation Coefficient ≤ 1.0 /// /// [flowArea] `cm²` /// /// [tamean] `cm/s` /// /// [coefficient] `[0.5,1.0]` /// /// return `cm³/s` @override double flowVolTAMAX(double flowArea, double tamean, double coefficient) { if (coefficient < 0.5 || coefficient > 1.0) { return double.nan; } double flowVol = flowArea * tamean.abs() * coefficient; return flowVol; } @override double pi(double ps, double mdEd, double tamax) { double pi = (ps - mdEd) / tamax; return pi.abs(); } @override double medianVelocity(double ps, double ed) { ps = ps.abs(); ed = ed.abs(); return (ps + 2 * ed) / 3; } @override double perimeter(double d1, double d2) { double perimeter = 0.0; if (!almostEquals(d1, 0.0) && !almostEquals(d2, 0.0)) { if (almostEquals(d1, d2)) { // Circle perimeter = d1 * math.pi; } else if (d1 < 0.001 * d2) { // Straight line perimeter = 2.0 * d2; } else if (d2 < 0.001 * d1) { // Straight line perimeter = 2.0 * d1; } else { // Ellipse double a = d1 > d2 ? d1 / 2.0 : d2 / 2.0; // longest half axis double b = d1 < d2 ? d1 / 2.0 : d2 / 2.0; // shortest half axis if (a < 0.0001) { // Long axis too short to make computations. Return zero. perimeter = 0.0; } else { // Calculate k, kc double k = math.sqrt(a * a - b * b) / a; double kc2 = 1 - k * k; perimeter = 4.0 * a * ellipseCircum(math.sqrt(kc2), kc2); } } } return perimeter; } bool almostEquals(double a, double b, [double epsilon = 1e-10]) { return (a - b).abs() < epsilon; } double ellipseCircum(double k, double kc2) { // Placeholder for the EllipseCircum method equivalent. // You can implement the complete elliptic integral of the second kind here // or use an approximation method for ellipse circumference. // For the sake of this example, we will use Ramanujan's approximation. double h = math.pow((k - 1), 2) / math.pow((k + 1), 2); return math.pi * (1 + 3 * h / (10 + math.sqrt(4 - 3 * h))); } @override double volume(double d) { double volume = 0.0; if (!almostEquals(d, 0.0)) { volume = math.pow(d, 3) * math.pi / 6.0; } return volume; } } class AnimalsGeneralFormulas extends BaseGeneralFormulas {}