klhztrader/KLHZ.Trader.Core.Math/Declisions/Utils/FFT.cs

using KLHZ.Trader.Core.Math.Declisions.Dtos.FFT;
using KLHZ.Trader.Core.Math.Declisions.Dtos.FFT.Enums;
using MathNet.Numerics;
using MathNet.Numerics.IntegralTransforms;
using static System.Runtime.InteropServices.JavaScript.JSType;

namespace KLHZ.Trader.Core.Math.Declisions.Utils
{
    public static class FFT
    {
        public static ValueAmplitudePosition Check(FFTAnalyzeResult fftData, DateTime timestamp)
        {
            var value = (decimal)CalcAmplitude(fftData.Harmonics, fftData.StartTime, timestamp);
            var value2 = (decimal)CalcExtremum(fftData.Harmonics, fftData.StartTime, timestamp);
            if (value > fftData.Upper30Decil)
            {
                return ValueAmplitudePosition.UpperThen30Decil;
            }
            else if (value < fftData.Mediana && System.Math.Sign(value2) >= 0)
            {
                return ValueAmplitudePosition.LowerThenMediana;
            }
            else
            {
                return ValueAmplitudePosition.Middle;
            }
        }

        public static FFTAnalyzeResult Analyze(DateTime[] timestamps, decimal[] values, string key, TimeSpan minPeriod, TimeSpan maxPeriod)
        {
            var harmonics = GetHarmonics(values, timestamps[timestamps.Length - 1] - timestamps[0], minPeriod, maxPeriod);
            var newValues = new decimal[timestamps.Length];
            var newValues2 = new decimal[timestamps.Length];
            var startTime = timestamps[0];
            for (int i = 0; i < timestamps.Length; i++)
            {
                newValues[i] = (decimal)CalcAmplitude(harmonics, startTime, timestamps[i]);
                newValues2[i] = (decimal)CalcExtremum(harmonics, startTime, timestamps[i]);
            }

            newValues = newValues.Order().ToArray();
            var ma = newValues2.Max();
            var mi = newValues2.Min();
            return new FFTAnalyzeResult()
            {
                Key = key,
                Harmonics = harmonics,
                LastTime = timestamps[timestamps.Length - 1],
                StartTime = startTime,
                Mediana = newValues[newValues.Length / 2],
                Upper30Decil = newValues[(int)(newValues.Length * 0.7)],
                Lower20Decil = newValues[(int)(newValues.Length * 0.2)],
                Max = newValues.Max(),
                Min = newValues.Min(),
                Length = values.Length,
            };
        }

        public static FFTAnalyzeResult ReAnalyze(FFTAnalyzeResult result, string key, TimeSpan minPeriod, TimeSpan maxPeriod)
        {
            var tmp = new List<Harmonic>();
            for (int i = 0; i < result.Harmonics.Length; i++)
            {
                var per = CaclHarmonycPeriod(result.LastTime - result.StartTime, result.Length, i);
                if (per >= minPeriod && per <= maxPeriod)
                {
                    tmp.Add(result.Harmonics[i]);
                }
            }

            var harms = tmp.ToArray();
            var newValues = new decimal[result.Length];
            var newValues2 = new decimal[result.Length];
            var time = result.StartTime;
            var dt = (result.LastTime - result.StartTime).TotalSeconds / result.Length;
            for (int i = 0; i < result.Length; i++)
            {
                var currentTime = time.AddSeconds(i* dt);
                newValues[i] = (decimal)CalcAmplitude(harms, result.StartTime, currentTime);
                newValues2[i] = (decimal)CalcExtremum(harms, result.StartTime, currentTime);
            }

            newValues = newValues.Order().ToArray();
            var ma = newValues2.Max();
            var mi = newValues2.Min();
            return new FFTAnalyzeResult()
            {
                Key = key,
                Harmonics = harms,
                LastTime = result.LastTime,
                StartTime = result.StartTime,
                Mediana = newValues[newValues.Length / 2],
                Upper30Decil = newValues[(int)(newValues.Length * 0.7)],
                Lower20Decil = newValues[(int)(newValues.Length * 0.2)],
                Max = newValues.Max(),
                Min = newValues.Min(),
            };
        }

        public static Harmonic[] GetHarmonics(decimal[] data, TimeSpan period, TimeSpan minPeriod, TimeSpan maxPeriod)
        {
            var arrv = data.Select(d => new Complex32((float)d, 0)).ToArray();
            Fourier.Forward(arrv);
            var res = new List<Harmonic>();
            for (int i = 1; i < arrv.Length; i++)
            {
                var per = CaclHarmonycPeriod(period, data.Length, i);
                if (per >= minPeriod && per <= maxPeriod)
                {
                    res.Add(new Harmonic()
                    {
                        Imaginary = arrv[i].Imaginary,
                        Real = arrv[i].Real,
                        Magnitude = arrv[i].Magnitude,
                        Period = per,
                        Phase = arrv[i].Phase,
                    });
                }
            }

            return res.ToArray();
        }

        public static double CalcAmplitude(Harmonic[] harmonics, DateTime startTime, DateTime currentTime)
        {
            var sum = 0d;
            var timeSpan = currentTime - startTime;
            for (int i = 0; i < harmonics.Length; i++)
            {
                var value = harmonics[i].Magnitude * System.Math.Cos(2 * System.Math.PI * timeSpan.TotalSeconds / harmonics[i].Period.TotalSeconds + harmonics[i].Phase);
                sum += value;
            }
            return sum;
        }

        public static double CalcExtremum(Harmonic[] harmonics, DateTime startTime, DateTime currentTime)
        {
            var sum = 0d;
            var timeSpan = currentTime - startTime;
            for (int i = 0; i < harmonics.Length; i++)
            {
                var value = -harmonics[i].Magnitude * System.Math.Sin(2 * System.Math.PI * timeSpan.TotalSeconds / harmonics[i].Period.TotalSeconds + harmonics[i].Phase);
                sum += value;
            }
            return sum;
        }

        internal static TimeSpan GetMainHarmonictPeriod(double[] data, TimeSpan period)
        {
            var arrv = data.Select(d => new Complex32((float)d, 0)).ToArray();
            Fourier.Forward(arrv);
            var res = new List<Harmonic>();
            var times = new TimeSpan[arrv.Length];
            var mags = arrv.Select(a => a.Magnitude).ToArray();
            var phases = arrv.Select(a => a.Phase).ToArray();
            var max = 0f;
            var kmax = 1;
            var sum = 0f;
            for (int i = 0; i < arrv.Length / 2; i++)
            {
                if (i == 0)
                {
                    res.Add(new Harmonic()
                    {
                        Magnitude = arrv[i].Magnitude,
                        Period = TimeSpan.MaxValue,
                        Phase = arrv[i].Phase,
                    });
                }
                else
                {
                    times[i] = CaclHarmonycPeriod(period, data.Length, i);
                    res.Add(new Harmonic()
                    {
                        Magnitude = arrv[i].Magnitude,
                        Period = times[i],
                        Phase = arrv[i].Phase,
                    });
                    var mag = arrv[i].Magnitude;
                    sum += mag;
                    if (mag > max)
                    {
                        max = mag;
                        kmax = i;
                    }
                }
            }
            return CaclHarmonycPeriod(period, data.Length, kmax);
        }

        internal static TimeSpan CaclHarmonycPeriod(TimeSpan signalLength, int signalLengthItems, int harmonyNumber)
        {
            var fdiscretisation = signalLengthItems / signalLength.TotalSeconds;
            var fharm = harmonyNumber * fdiscretisation / signalLengthItems;
            return TimeSpan.FromSeconds(1 / fharm);
        }

        internal static double CalcCurrentPhase(Complex32[] spectr, int harmonyNumber)
        {
            var item = spectr[harmonyNumber];
            return System.Math.Atan(item.Imaginary / item.Real);
        }

        internal static (double min, double max) CalcPhaseRangeFoxMax(double aSin, double aCos, double initPhase, double level)
        {
            return CalcPhaseRange(aSin, aCos, initPhase, level, CheckMaxValue);
        }

        internal static (double min, double max) CalcPhaseRangeFoxMin(double aSin, double aCos, double initPhase, double level)
        {
            return CalcPhaseRange(aSin, aCos, initPhase, level, CheckMinValue);
        }

        internal static (double min, double max) CalcPhaseRange(double aSin, double aCos, double initPhase, double level, Func<double, double, double, double, bool> comparer)
        {
            var x = new List<double>();
            var xIndexes = new List<int>();
            var y = new List<double>();
            var max = double.MinValue;
            var min = double.MaxValue;
            var _2pi = 2 * System.Math.PI;
            for (double i = 0; i <= 10 * System.Math.PI; i += 0.01)
            {
                var df = (((i) / _2pi) % 1) * _2pi;
                var val = aSin * System.Math.Sin(i + initPhase) + aCos * System.Math.Cos(i + initPhase);
                if (val > max)
                {
                    max = val;
                }
                if (val < min)
                {
                    min = val;
                }
                x.Add(df);
                y.Add(val);
            }

            int start = -2;
            int prevI = -2;
            int end = -2;

            var drange = -2;

            var minPhaseTmp = 0d;
            var maxPhaseTmp = 0d;

            var minPhase = 0d;
            var maxPhase = 0d;
            for (int i = 0; i < y.Count; i++)
            {
                if (comparer(y[i], min, max, level))
                {
                    if (start < 0)
                    {
                        minPhaseTmp = x[i];
                        start = i;
                    }
                }
                else if (start >= 0 && i - prevI == 1)
                {
                    end = prevI;
                    maxPhaseTmp = x[end];
                    var drangeTmp = end - start;
                    if (drangeTmp > drange)
                    {
                        drange = drangeTmp;
                        minPhase = minPhaseTmp;
                        maxPhase = maxPhaseTmp;
                    }
                    start = -2;
                }
                prevI = i;
            }


            return (minPhase, maxPhase);
        }

        internal static bool CheckMaxValue(double val, double min, double max, double relativeValue)
        {
            return (val > (1 - relativeValue) * max);
        }

        internal static bool CheckMinValue(double val, double min, double max, double relativeValue)
        {
            return (val < (1 - relativeValue) * min);
        }
    }
}