ningshuxia
2024-04-28 cc3788df309c28a19f61331a5ec5379717799a9b
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
namespace IStation.Curve
{
    /// <summary>
    /// 样条曲线(四次拟合):使用MATH.net
    /// </summary>
    public class TwinRelateCurve4M : IFitPoint
    {
        /// <summary>
        /// 
        /// </summary>
        public TwinRelateCurve4M(CurveExpress express)
        {
            _express = express;
        }
 
        /// <summary>
        /// 
        /// </summary>
        public TwinRelateCurve4M(List<CurvePoint> points)
        {
            if (points == null || points.Count < 5)
            {
                throw new Exception("出现数值计算问题 ,无法解出四次样条曲线!");
            }
 
            _express = BuildCurveExpress(points);
        }
 
        //曲线表达式
        private CurveExpress _express = null;
 
 
        /// <summary>
        /// 获取拟合点Y
        /// </summary>
        public double GetFitPointY(double x)
        {
            if (_express == null)
                return default;
            //不用减掉最小值, 和alglib有点不一样
            return _express.Index4 * x * x * x * x + _express.Index3 * x * x * x + _express.Index2 * x * x + _express.Index1 * x + _express.Index0;
        }
 
        /// <summary>
        /// 获取拟合点列表
        /// </summary>
        public List<CurvePoint> GetFitPoints(int pointNumber)
        {
            if (_express == null)
                return default;
            return GetFitPointsByXRange(_express.Min, _express.Max, pointNumber);
        }
 
        /// <summary>
        /// 通过区间获取拟合点列表
        /// </summary>
        public List<CurvePoint> GetFitPointsByXRange(double minX, double maxX, int pointNumber)
        {
            if (_express == null)
                return default;
 
            if (pointNumber < 1)
                pointNumber = 12;
            if (pointNumber > 10000)
                pointNumber = 10000;
            double space = (maxX - minX) / (pointNumber - 1);
            if (space < 0.0001)
                return default;
            var points = new List<CurvePoint>();
            for (int i = 0; i < pointNumber; i++)
            {
                double x = space * i + minX;
                double y = this.GetFitPointY(x);
                points.Add(new CurvePoint(x, y));
            }
            return points;
        }
 
        /// <summary>
        /// 创建曲线表达式
        /// </summary>
        public static CurveExpress BuildCurveExpress(List<CurvePoint> points)
        {
            if (points == null || points.Count < 5)
                return default;
 
            points = points.OrderBy(x => x.X).ToList();
            CurveExpress express = new CurveExpress();
            express.FitType = eFitType.FourM;
            express.DefinePoints = points.Select(x => new CurvePoint(x)).ToList();
            express.Min = points.First().X;
            express.Max = points.Last().X;
 
            var X = (from x in points select x.X).ToArray();
            var Y = (from x in points select x.Y).ToArray();
            var res = MathNet.Numerics.Fit.Polynomial(X, Y, 4);
 
            //0次方
            express.Index0 = res[0];
            //1次方
            express.Index1 = res[1];
            //2次方
            express.Index2 = res[2];
            //3次方
            express.Index3 = res[3];
            //4次方
            express.Index4 = res[4];
 
            return express;
        }
 
    }
}