From 967096d4f5084280465d37352da40fdd52af79dd Mon Sep 17 00:00:00 2001 From: Miha Kralj Date: Thu, 4 Dec 2025 13:49:05 -0800 Subject: [PATCH] Refactor and optimize various components of QuanTAlib - Removed WmaVector class to streamline weighted moving average calculations. - Simplified RingBuffer implementation by removing unnecessary comments and improving clarity. - Enhanced SIMD extensions for better performance and readability. - Updated TBar and TBarSeries classes to improve property calculations and reduce overhead. - Cleaned up TValue struct by removing redundant comments. - Added comprehensive unit tests for IndicatorExtensions and TrimaIndicator to ensure functionality and correctness. --- lib/averages/ema/Ema.Notebook.dib | 81 +---- lib/averages/ema/Ema.cs | 68 +--- lib/averages/ema/Ema.md | 28 +- lib/averages/ema/EmaVector.Tests.cs | 378 ---------------------- lib/averages/ema/EmaVector.cs | 327 ------------------- lib/averages/sma/Sma.Notebook.dib | 84 +---- lib/averages/sma/Sma.cs | 100 +----- lib/averages/sma/Sma.md | 28 +- lib/averages/sma/SmaVector.Tests.cs | 371 --------------------- lib/averages/sma/SmaVector.cs | 185 ----------- lib/averages/trima/Trima.Tests.cs | 5 +- lib/averages/trima/Trima.cs | 126 ++------ lib/averages/trima/TrimaVector.Tests.cs | 363 --------------------- lib/averages/trima/TrimaVector.cs | 241 -------------- lib/averages/wma/Wma.Notebook.dib | 64 +--- lib/averages/wma/Wma.cs | 235 ++------------ lib/averages/wma/Wma.md | 28 +- lib/averages/wma/WmaVector.Tests.cs | 407 ------------------------ lib/averages/wma/WmaVector.cs | 264 --------------- lib/core/ringbuffer/RingBuffer.cs | 10 +- lib/core/simd/SimdExtensions.cs | 31 +- lib/core/tbar/tbar.cs | 4 +- lib/core/tbarseries/tbarseries.cs | 7 - lib/core/tseries/tseries.cs | 18 +- lib/core/tvalue/tvalue.cs | 10 - quantower/IndicatorExtensions.Tests.cs | 148 +++++++++ quantower/TrimaIndicator.Tests.cs | 143 +++++++++ 27 files changed, 387 insertions(+), 3367 deletions(-) delete mode 100644 lib/averages/ema/EmaVector.Tests.cs delete mode 100644 lib/averages/ema/EmaVector.cs delete mode 100644 lib/averages/sma/SmaVector.Tests.cs delete mode 100644 lib/averages/sma/SmaVector.cs delete mode 100644 lib/averages/trima/TrimaVector.Tests.cs delete mode 100644 lib/averages/trima/TrimaVector.cs delete mode 100644 lib/averages/wma/WmaVector.Tests.cs delete mode 100644 lib/averages/wma/WmaVector.cs create mode 100644 quantower/IndicatorExtensions.Tests.cs create mode 100644 quantower/TrimaIndicator.Tests.cs diff --git a/lib/averages/ema/Ema.Notebook.dib b/lib/averages/ema/Ema.Notebook.dib index bcf534e1..4818ac63 100644 --- a/lib/averages/ema/Ema.Notebook.dib +++ b/lib/averages/ema/Ema.Notebook.dib @@ -16,7 +16,7 @@ This notebook demonstrates: 1. **Manual Data Processing**: Understanding Batch vs. Streaming modes. 2. **Streaming with `isNew`**: Handling intra-bar updates. 3. **Large Dataset Processing**: Using Geometric Brownian Motion (GBM) generated data. -4. **Vectorized Operations**: Calculating multiple EMAs simultaneously. +4. **Handling Invalid Values**: Last-value substitution for NaN/Infinity. #!csharp @@ -165,63 +165,9 @@ Console.WriteLine($"Streaming Last Value: {lastStreamVal.Value:F2}"); #!markdown -## 4. Vectorized EMA (Multiple Periods) +## 4. Handling Invalid Values (NaN/Infinity) -`EmaVector` allows calculating multiple EMAs (e.g., 9, 12, 26) simultaneously. This is optimized for performance using SIMD where available. - -### Vectorized Batch - -#!csharp - -int[] periods = { 9, 12, 26 }; -Console.WriteLine($"\n--- Vectorized Batch EMA (Periods: {string.Join(", ", periods)}) ---"); - -var emaVectorBatch = new EmaVector(periods); -var vectorBatchResults = emaVectorBatch.Calculate(closeSeries); - -for (int i = 0; i < periods.Length; i++) -{ - Console.WriteLine($"EMA({periods[i]}) Last Value: {vectorBatchResults[i].Last().Value:F2}"); -} - -#!markdown - -### Vectorized Streaming - -#!csharp - -Console.WriteLine($"\n--- Vectorized Streaming EMA (Periods: {string.Join(", ", periods)}) ---"); - -var emaVectorStream = new EmaVector(periods); -TValue[] lastVectorVal = null; - -foreach(var item in closeSeries) -{ - lastVectorVal = emaVectorStream.Update(item); -} - -for (int i = 0; i < periods.Length; i++) -{ - Console.WriteLine($"EMA({periods[i]}) Last Value: {lastVectorVal[i].Value:F2}"); -} - -// Verification -bool allMatch = true; -for (int i = 0; i < periods.Length; i++) -{ - if (Math.Abs(vectorBatchResults[i].Last().Value - lastVectorVal[i].Value) > 1e-10) - { - allMatch = false; - break; - } -} -Console.WriteLine($"\nAll Vectorized Stream/Batch values match: {allMatch}"); - -#!markdown - -## 5. Handling Invalid Values (NaN/Infinity) - -Both `Ema` and `EmaVector` use **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. +`Ema` uses **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. #!csharp @@ -271,24 +217,3 @@ for (int i = 0; i < seriesWithNaN.Count; i++) var inputStr = double.IsFinite(input) ? input.ToString("F2") : input.ToString(); Console.WriteLine($" {inputStr,-10} → {output:F2} (IsFinite: {double.IsFinite(output)})"); } - -#!csharp - -Console.WriteLine("\n--- Vectorized EMA with Invalid Values ---"); - -int[] periodsNaN = { 5, 10 }; -var emaVectorNaN = new EmaVector(periodsNaN); - -// Feed values including invalid ones -var inputsNaN = new double[] { 100, 110, double.NaN, 120, double.PositiveInfinity, 130 }; -var time = DateTime.Now; - -foreach (var val in inputsNaN) -{ - var results = emaVectorNaN.Update(new TValue(time, val)); - var inputStr = double.IsFinite(val) ? val.ToString("F2") : val.ToString(); - Console.WriteLine($"Input: {inputStr,-10} → EMA(5): {results[0].Value:F2}, EMA(10): {results[1].Value:F2}"); - time = time.AddMinutes(1); -} - -Console.WriteLine("\nAll outputs are finite - invalid inputs were substituted with last valid values."); diff --git a/lib/averages/ema/Ema.cs b/lib/averages/ema/Ema.cs index d1a539fa..c843d073 100644 --- a/lib/averages/ema/Ema.cs +++ b/lib/averages/ema/Ema.cs @@ -25,30 +25,18 @@ namespace QuanTAlib; /// public class Ema { - private struct State : IEquatable + private struct State { public double Ema; - public double E; // Compensator: decays from 1.0 to 1e-10 for bias correction - public bool IsHot; // True when 95% coverage reached (E <= 0.05) - public bool IsCompensated; // True when compensator fully decayed (E <= 1e-10) + public double E; + public bool IsHot; + public bool IsCompensated; public static State New() => new() { Ema = 0, E = 1.0, IsHot = false, IsCompensated = false }; - - public readonly bool Equals(State other) => - Ema == other.Ema && E == other.E && IsHot == other.IsHot && IsCompensated == other.IsCompensated; - - public override readonly bool Equals(object? obj) => - obj is State other && Equals(other); - - public override readonly int GetHashCode() => - HashCode.Combine(Ema, E, IsHot, IsCompensated); - - public static bool operator ==(State left, State right) => left.Equals(right); - public static bool operator !=(State left, State right) => !left.Equals(right); } private readonly double _alpha; - private readonly double _decay; // Pre-calculated (1.0 - alpha) to avoid subtraction per tick + private readonly double _decay; private State _state = State.New(); private State _p_state = State.New(); private double _lastValidValue; @@ -97,9 +85,6 @@ public class Ema /// public bool IsHot => _state.IsHot; - /// - /// Gets a valid input value, using last-value substitution for non-finite inputs. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private double GetValidValue(double input) { @@ -111,17 +96,9 @@ public class Ema return _lastValidValue; } - // 95% coverage threshold: E = 1 - 0.95 = 0.05 private const double COVERAGE_THRESHOLD = 0.05; - // Compensator decay threshold for bias correction private const double COMPENSATOR_THRESHOLD = 1e-10; - /// - /// Core EMA calculation kernel. - /// Assumes input has already been validated via GetValidValue(). - /// IsHot becomes true at 95% coverage (E <= 0.05). - /// Bias correction continues until compensator decays to 1e-10. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private static double Compute(double input, double alpha, double decay, ref State state) { @@ -132,11 +109,9 @@ public class Ema { state.E *= decay; - // IsHot triggers at 95% coverage if (!state.IsHot && state.E <= COVERAGE_THRESHOLD) state.IsHot = true; - // Continue bias correction until compensator fully decays if (state.E <= COMPENSATOR_THRESHOLD) { state.IsCompensated = true; @@ -155,9 +130,6 @@ public class Ema return result; } - /// - /// Core calculation kernel that handles both batch and streaming-continuation. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private static void CalculateCore(ReadOnlySpan source, Span output, double alpha, ref State state, ref double lastValidValue) { @@ -165,8 +137,6 @@ public class Ema double decay = 1.0 - alpha; int i = 0; - // Phase 1: Warmup with bias correction - // If state is already compensated, this loop is skipped if (!state.IsCompensated) { for (; i < len && state.E > COMPENSATOR_THRESHOLD; i++) @@ -189,7 +159,6 @@ public class Ema state.IsCompensated = true; } - // Phase 2: Hot loop for (; i < len; i++) { double val = source[i]; @@ -203,12 +172,6 @@ public class Ema } } - /// - /// Updates EMA with the given value. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Compensated EMA value [MethodImpl(MethodImplOptions.AggressiveInlining)] public TValue Update(TValue input, bool isNew = true) { @@ -221,19 +184,12 @@ public class Ema _state = _p_state; } - // Last-value substitution: replace non-finite inputs with last valid value double val = GetValidValue(input.Value); val = Compute(val, _alpha, _decay, ref _state); Value = new TValue(input.Time, val); return Value; } - /// - /// Updates EMA with the entire series. - /// Uses split-loop optimization: warmup phase with bias correction, then branchless hot loop. - /// - /// Input series - /// EMA series public TSeries Update(TSeries source) { if (source.Count == 0) return new TSeries(new List(), new List()); @@ -249,9 +205,6 @@ public class Ema var sourceValues = source.Values; var sourceTimes = source.Times; - // 1. Fast Batch Calculation - // Uses the unified CalculateCore to handle both new and continuing states - // Optimization: Copy state to locals to allow JIT register allocation State state = _state; double lastValidValue = _lastValidValue; @@ -260,7 +213,6 @@ public class Ema _state = state; _lastValidValue = lastValidValue; - // Copy Times sourceTimes.CopyTo(tSpan); _p_state = _state; @@ -299,15 +251,6 @@ public class Ema Calculate(source, output, alpha); } - /// - /// Calculates EMA in-place using alpha, writing results to pre-allocated output span. - /// Zero-allocation method for maximum performance. - /// Bias correction continues until compensator decays to 1e-10. - /// Uses split-loop optimization: warmup phase with bias correction, then branchless hot loop. - /// - /// Input values - /// Output span (must be same length as source) - /// Smoothing factor (0 < alpha <= 1) [MethodImpl(MethodImplOptions.AggressiveInlining)] public static void Calculate(ReadOnlySpan source, Span output, double alpha) { @@ -318,7 +261,6 @@ public class Ema if (source.Length == 0) return; - // Initialize default state for static calculation State state = State.New(); double lastValid = 0; diff --git a/lib/averages/ema/Ema.md b/lib/averages/ema/Ema.md index 1cd287eb..c45ad1fe 100644 --- a/lib/averages/ema/Ema.md +++ b/lib/averages/ema/Ema.md @@ -115,33 +115,9 @@ Console.WriteLine($"Last EMA: {emaOutput[^1]}"); * **Hunter's bias correction**: Same accuracy as TSeries API * **Compatible** with `ArrayPool` for buffer management -### Multi-Alpha EMA (`EmaVector`) - -The `EmaVector` class is a SIMD-optimized implementation for calculating multiple EMAs with different periods on the same input series simultaneously. It leverages hardware intrinsics (AVX/SSE) for high performance. - -```csharp -using QuanTAlib; - -// Initialize with multiple periods -int[] periods = { 9, 12, 26 }; -var emaVector = new EmaVector(periods); - -// Streaming update -TValue[] results = emaVector.Update(new TValue(time, price)); - -// Access values -Console.WriteLine($"EMA(9): {results[0].Value}"); -Console.WriteLine($"EMA(12): {results[1].Value}"); -Console.WriteLine($"EMA(26): {results[2].Value}"); - -// Batch calculation -TSeries source = ...; -TSeries[] seriesResults = emaVector.Calculate(source); -``` - ### Handling Invalid Values (NaN/Infinity) -Both `Ema` and `EmaVector` use **last-value substitution** for handling invalid inputs: +`Ema` uses **last-value substitution** for handling invalid inputs: ```csharp var ema = new Ema(10); @@ -166,13 +142,11 @@ var results = ema.Update(series); // All values are finite * When `NaN`, `PositiveInfinity`, or `NegativeInfinity` is encountered, the last valid value is substituted * This provides output continuity instead of propagating invalid values -* Both scalar (`Ema`) and SIMD (`EmaVector`) implementations use identical logic * `Reset()` clears the last valid value, so the next valid input establishes a new baseline ### Performance Characteristics * **O(1) Complexity:** The calculation time is constant regardless of the period length. -* **SIMD Optimization:** `EmaVector` processes multiple periods in parallel using vector instructions, significantly reducing CPU cycles for multi-timeframe analysis. * **Zero Allocation:** The streaming `Update` method is designed to be allocation-free (excluding the return struct). ## Interpretation Details diff --git a/lib/averages/ema/EmaVector.Tests.cs b/lib/averages/ema/EmaVector.Tests.cs deleted file mode 100644 index ab466839..00000000 --- a/lib/averages/ema/EmaVector.Tests.cs +++ /dev/null @@ -1,378 +0,0 @@ - -namespace QuanTAlib.Tests; - -public class EmaVectorTests -{ - [Fact] - public void Initialization_WithPeriods_SetsCorrectAlphas() - { - int[] periods = { 10, 20 }; - var emaVector = new EmaVector(periods); - - var res = emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(100.0, res[0].Value, 1e-9); - Assert.Equal(100.0, res[1].Value, 1e-9); - } - - [Fact] - public void Initialization_WithAlphas_Works() - { - double[] alphas = { 0.1, 0.2, 0.5 }; - var emaVector = new EmaVector(alphas); - - var res = emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(3, res.Length); - Assert.Equal(100.0, res[0].Value, 1e-9); - Assert.Equal(100.0, res[1].Value, 1e-9); - Assert.Equal(100.0, res[2].Value, 1e-9); - } - - [Fact] - public void Initialization_WithZeroPeriod_ThrowsArgumentException() - { - int[] periods = { 10, 0, 20 }; - - Assert.Throws(() => new EmaVector(periods)); - } - - [Fact] - public void Initialization_WithNegativePeriod_ThrowsArgumentException() - { - int[] periods = { 10, -5, 20 }; - - Assert.Throws(() => new EmaVector(periods)); - } - - [Fact] - public void Initialization_WithZeroAlpha_ThrowsArgumentException() - { - double[] alphas = { 0.1, 0.0, 0.5 }; - - Assert.Throws(() => new EmaVector(alphas)); - } - - [Fact] - public void Initialization_WithNegativeAlpha_ThrowsArgumentException() - { - double[] alphas = { 0.1, -0.1, 0.5 }; - - Assert.Throws(() => new EmaVector(alphas)); - } - - [Fact] - public void Initialization_WithAlphaGreaterThanOne_ThrowsArgumentException() - { - double[] alphas = { 0.1, 1.5, 0.5 }; - - Assert.Throws(() => new EmaVector(alphas)); - } - - [Fact] - public void Initialization_WithAlphaEqualToOne_Works() - { - double[] alphas = { 0.1, 1.0, 0.5 }; - var emaVector = new EmaVector(alphas); - - var res = emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(3, res.Length); - } - - [Fact] - public void Calc_Streaming_MatchesSingleEma() - { - int[] periods = { 5, 10, 20 }; - var emaVector = new EmaVector(periods); - var emaSingles = periods.Select(p => new Ema(p)).ToArray(); - - var values = new double[] { 10, 20, 30, 40, 50, 40, 30, 20, 10 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = emaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = emaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - Assert.Equal(singleRes.Time, multiRes[i].Time); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Calc_Series_MatchesSingleEma() - { - int[] periods = { 5, 10, 20 }; - var emaVector = new EmaVector(periods); - var emaSingles = periods.Select(p => new Ema(p)).ToArray(); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var multiRes = emaVector.Calculate(series); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = emaSingles[i].Update(series); - - Assert.Equal(singleRes.Count, multiRes[i].Count); - for (int j = 0; j < len; j++) - { - Assert.Equal(singleRes.Values[j], multiRes[i].Values[j], 1e-8); - } - } - } - - [Fact] - public void Calc_Series_MatchesStreaming() - { - int[] periods = { 5, 10, 20 }; - var emaVectorBatch = new EmaVector(periods); - var emaVectorStream = new EmaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var batchRes = emaVectorBatch.Calculate(series); - - for (int i = 0; i < len; i++) - { - var tVal = new TValue(new DateTime(t[i], DateTimeKind.Utc), v[i]); - var streamRes = emaVectorStream.Update(tVal); - - for (int j = 0; j < periods.Length; j++) - { - Assert.Equal(batchRes[j].Values[i], streamRes[j].Value, 1e-9); - } - } - } - - [Fact] - public void Calculate_Static_MatchesInstanceMethod() - { - int[] periods = { 5, 10, 20 }; - - int len = 50; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var instanceEma = new EmaVector(periods); - var instanceRes = instanceEma.Calculate(series); - - var staticRes = EmaVector.Calculate(series, periods); - - for (int i = 0; i < periods.Length; i++) - { - Assert.Equal(instanceRes[i].Count, staticRes[i].Count); - for (int j = 0; j < len; j++) - { - Assert.Equal(instanceRes[i].Values[j], staticRes[i].Values[j], 1e-9); - } - } - } - - [Fact] - public void Reset_ClearsState() - { - int[] periods = { 10 }; - var emaVector = new EmaVector(periods); - - emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - emaVector.Reset(); - - var res = emaVector.Update(new TValue(DateTime.UtcNow, 200.0)); - - Assert.Equal(200.0, res[0].Value, 1e-9); - } - - [Fact] - public void Update_NaN_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var emaVector = new EmaVector(periods); - - emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - emaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterNaN = emaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in resultAfterNaN) - { - Assert.True(double.IsFinite(result.Value), $"Expected finite value but got {result.Value}"); - } - } - - [Fact] - public void Update_Infinity_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var emaVector = new EmaVector(periods); - - emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - emaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterPosInf = emaVector.Update(new TValue(DateTime.UtcNow, double.PositiveInfinity)); - foreach (var result in resultAfterPosInf) - { - Assert.True(double.IsFinite(result.Value)); - } - - var resultAfterNegInf = emaVector.Update(new TValue(DateTime.UtcNow, double.NegativeInfinity)); - foreach (var result in resultAfterNegInf) - { - Assert.True(double.IsFinite(result.Value)); - } - } - - [Fact] - public void Update_MultipleNaN_ContinuesWithLastValid() - { - int[] periods = { 5, 10 }; - var emaVector = new EmaVector(periods); - - emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - emaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - emaVector.Update(new TValue(DateTime.UtcNow, 120.0)); - - var r1 = emaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r2 = emaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r3 = emaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in r1) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r2) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r3) Assert.True(double.IsFinite(result.Value)); - } - - [Fact] - public void Calculate_Series_HandlesNaN() - { - int[] periods = { 5, 10 }; - var emaVector = new EmaVector(periods); - - var t = new System.Collections.Generic.List(); - var v = new System.Collections.Generic.List(); - var now = DateTime.UtcNow; - - t.Add(now.Ticks); v.Add(100.0); - t.Add(now.AddMinutes(1).Ticks); v.Add(110.0); - t.Add(now.AddMinutes(2).Ticks); v.Add(double.NaN); - t.Add(now.AddMinutes(3).Ticks); v.Add(120.0); - t.Add(now.AddMinutes(4).Ticks); v.Add(double.PositiveInfinity); - t.Add(now.AddMinutes(5).Ticks); v.Add(130.0); - - var series = new TSeries(t, v); - var results = emaVector.Calculate(series); - - foreach (var periodResults in results) - { - foreach (var val in periodResults.Values) - { - Assert.True(double.IsFinite(val), $"Expected finite value but got {val}"); - } - } - } - - [Fact] - public void Reset_ClearsLastValidValue() - { - int[] periods = { 10 }; - var emaVector = new EmaVector(periods); - - emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - emaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - emaVector.Reset(); - - var result = emaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - Assert.Equal(50.0, result[0].Value, 1e-9); - } - - [Fact] - public void NaN_Handling_MatchesSingleEma() - { - int[] periods = { 5, 10, 20 }; - var emaVector = new EmaVector(periods); - var emaSingles = periods.Select(p => new Ema(p)).ToArray(); - - var values = new double[] { 10, 20, double.NaN, 40, double.PositiveInfinity, 60, 70 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = emaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = emaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Values_Property_UpdatesAfterUpdate() - { - int[] periods = { 5, 10 }; - var emaVector = new EmaVector(periods); - - var result = emaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(result[0].Value, emaVector.Values[0].Value); - Assert.Equal(result[1].Value, emaVector.Values[1].Value); - } - - [Fact] - public void Values_Property_UpdatesAfterCalculate() - { - int[] periods = { 5, 10 }; - var emaVector = new EmaVector(periods); - - var t = new System.Collections.Generic.List { 100, 200, 300 }; - var v = new System.Collections.Generic.List { 10.0, 20.0, 30.0 }; - var series = new TSeries(t, v); - - var results = emaVector.Calculate(series); - - Assert.Equal(results[0].Last.Value, emaVector.Values[0].Value, 1e-9); - Assert.Equal(results[1].Last.Value, emaVector.Values[1].Value, 1e-9); - } -} diff --git a/lib/averages/ema/EmaVector.cs b/lib/averages/ema/EmaVector.cs deleted file mode 100644 index 6db798c4..00000000 --- a/lib/averages/ema/EmaVector.cs +++ /dev/null @@ -1,327 +0,0 @@ -using System.Numerics; -using System.Runtime.CompilerServices; -using System.Runtime.InteropServices; - -namespace QuanTAlib; - -/// -/// Multi-Alpha Exponential Moving Average (EMA) - SIMD optimized. -/// Calculates multiple EMAs with different periods/alphas for the same input series in parallel. -/// Uses last-value substitution for invalid inputs (NaN/Infinity). -/// -[SkipLocalsInit] -public class EmaVector -{ - private readonly double[] _alphas; - private readonly double[] _emas; - private readonly double[] _Es; - private readonly double[] _p_emas; - private readonly double[] _p_Es; - private readonly int _count; - private double _lastValidValue; - - /// - /// Current EMA values for all periods. - /// - public ReadOnlySpan Values => _values; - - private readonly TValue[] _values; - - /// - /// Initializes EmaVector with specified periods. - /// - /// Array of periods - public EmaVector(int[] periods) - { - _count = periods.Length; - _alphas = new double[_count]; - _emas = new double[_count]; - _Es = new double[_count]; - _p_emas = new double[_count]; - _p_Es = new double[_count]; - _values = new TValue[_count]; - - for (int i = 0; i < _count; i++) - { - ArgumentOutOfRangeException.ThrowIfLessThanOrEqual(periods[i], 0); - _alphas[i] = 2.0 / (periods[i] + 1); - ResetAt(i); - } - } - - /// - /// Initializes EmaVector with specified alphas. - /// - /// Array of alphas - public EmaVector(double[] alphas) - { - _count = alphas.Length; - _alphas = new double[_count]; - _emas = new double[_count]; - _Es = new double[_count]; - _p_emas = new double[_count]; - _p_Es = new double[_count]; - _values = new TValue[_count]; - - for (int i = 0; i < _count; i++) - { - if (alphas[i] <= 0 || alphas[i] > 1) - throw new ArgumentOutOfRangeException(nameof(alphas), alphas[i], "Alpha must be between 0 (exclusive) and 1 (inclusive)"); - _alphas[i] = alphas[i]; - ResetAt(i); - } - } - - private void ResetAt(int index) - { - _emas[index] = 0.0; - _Es[index] = 1.0; - } - - /// - /// Gets a valid input value, using last-value substitution for non-finite inputs. - /// - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private double GetValidValue(double input) - { - if (double.IsFinite(input)) - { - _lastValidValue = input; - return input; - } - return _lastValidValue; - } - - /// - /// Resets all EMA states. - /// - public void Reset() - { - for (int i = 0; i < _count; i++) - { - ResetAt(i); - } - _lastValidValue = 0; - Array.Clear(_values); - } - - /// - /// Updates EMAs with the given value. - /// Uses last-value substitution: invalid inputs (NaN/Infinity) are replaced with - /// the last known good value, providing continuity in the output series. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Array of compensated EMA values - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public TValue[] Update(TValue input, bool isNew = true) - { - if (isNew) - { - Array.Copy(_emas, _p_emas, _count); - Array.Copy(_Es, _p_Es, _count); - } - else - { - Array.Copy(_p_emas, _emas, _count); - Array.Copy(_p_Es, _Es, _count); - } - - // Last-value substitution: replace non-finite inputs with last valid value - double val = GetValidValue(input.Value); - - // SIMD Loop - int vecCount = Vector.Count; - int i = 0; - - if (Vector.IsHardwareAccelerated && _count >= vecCount) - { - var vecInput = new Vector(val); - var vecOne = Vector.One; - var vecEpsilon = new Vector(1e-10); - - for (; i <= _count - vecCount; i += vecCount) - { - // Load state - var vecAlpha = new Vector(_alphas, i); - var vecEma = new Vector(_emas, i); - var vecE = new Vector(_Es, i); - - // Update EMA - // ema += alpha * (input - ema) - vecEma += vecAlpha * (vecInput - vecEma); - - // Update E (warmup factor) - // E *= (1 - alpha) - vecE *= (vecOne - vecAlpha); - - // Calculate compensated result - // res = ema / (1 - E) - var vecCompensated = vecEma / (vecOne - vecE); - - // Check warmup condition: E <= 1e-10 means "hot" (use raw EMA) - // Vector.LessThanOrEqual returns Vector with all-1s for true, all-0s for false - // We reinterpret as Vector for use with ConditionalSelect - var isHotMask = Vector.LessThanOrEqual(vecE, vecEpsilon); - - // Select result: if hot (E <= epsilon), use raw EMA; otherwise use compensated - // ConditionalSelect: mask=true -> first arg, mask=false -> second arg - var vecResult = Vector.ConditionalSelect( - Vector.AsVectorDouble(isHotMask), - vecEma, // Hot: use raw EMA - vecCompensated // Cold: use compensated - ); - - // Store state - vecEma.CopyTo(_emas, i); - vecE.CopyTo(_Es, i); - - // Store result - for (int j = 0; j < vecCount; j++) - { - _values[i + j] = new TValue(input.Time, vecResult[j]); - } - } - } - - // Scalar fallback for remaining items - for (; i < _count; i++) - { - double alpha = _alphas[i]; - _emas[i] += alpha * (val - _emas[i]); - - double result = _emas[i]; - if (_Es[i] > 1e-10) - { - _Es[i] *= (1.0 - alpha); - if (_Es[i] > 1e-10) - { - result = _emas[i] / (1.0 - _Es[i]); - } - } - - _values[i] = new TValue(input.Time, result); - } - - return _values; - } - - /// - /// Calculates EMAs for the entire series. - /// - /// Input series - /// Array of EMA series - public TSeries[] Calculate(TSeries source) - { - int len = source.Count; - var resultSeries = new TSeries[_count]; - - // Pre-allocate lists - var tLists = new List[_count]; - var vLists = new List[_count]; - - for (int i = 0; i < _count; i++) - { - tLists[i] = new List(len); - vLists[i] = new List(len); - CollectionsMarshal.SetCount(tLists[i], len); - CollectionsMarshal.SetCount(vLists[i], len); - } - - var sourceValues = source.Values; - var sourceTimes = source.Times; - - int vecCount = Vector.Count; - var vecOne = Vector.One; - var vecEpsilon = new Vector(1e-10); - - for (int t = 0; t < len; t++) - { - double val = sourceValues[t]; - long time = sourceTimes[t]; - - // Last-value substitution: replace non-finite inputs with last valid value - val = GetValidValue(val); - - var vecInput = new Vector(val); - - int i = 0; - if (Vector.IsHardwareAccelerated && _count >= vecCount) - { - for (; i <= _count - vecCount; i += vecCount) - { - var vecAlpha = new Vector(_alphas, i); - var vecEma = new Vector(_emas, i); - var vecE = new Vector(_Es, i); - - vecEma += vecAlpha * (vecInput - vecEma); - vecE *= (vecOne - vecAlpha); - - var vecCompensated = vecEma / (vecOne - vecE); - - // Check warmup condition: E <= 1e-10 means "hot" (use raw EMA) - var isHotMask = Vector.LessThanOrEqual(vecE, vecEpsilon); - - // Select result: if hot, use raw EMA; otherwise use compensated - var vecResult = Vector.ConditionalSelect( - Vector.AsVectorDouble(isHotMask), - vecEma, // Hot: use raw EMA - vecCompensated // Cold: use compensated - ); - - vecEma.CopyTo(_emas, i); - vecE.CopyTo(_Es, i); - - // Scatter results to lists - for (int j = 0; j < vecCount; j++) - { - CollectionsMarshal.AsSpan(tLists[i + j])[t] = time; - CollectionsMarshal.AsSpan(vLists[i + j])[t] = vecResult[j]; - } - } - } - - for (; i < _count; i++) - { - double alpha = _alphas[i]; - _emas[i] += alpha * (val - _emas[i]); - - double result = _emas[i]; - if (_Es[i] > 1e-10) - { - _Es[i] *= (1.0 - alpha); - if (_Es[i] > 1e-10) - { - result = _emas[i] / (1.0 - _Es[i]); - } - } - - CollectionsMarshal.AsSpan(tLists[i])[t] = time; - CollectionsMarshal.AsSpan(vLists[i])[t] = result; - } - } - - // Create TSeries and update Values - for (int i = 0; i < _count; i++) - { - resultSeries[i] = new TSeries(tLists[i], vLists[i]); - var lastT = CollectionsMarshal.AsSpan(tLists[i])[len - 1]; - var lastV = CollectionsMarshal.AsSpan(vLists[i])[len - 1]; - _values[i] = new TValue(lastT, lastV); - } - - return resultSeries; - } - - /// - /// Calculates EMAs for the entire series using specified periods. - /// - /// Input series - /// Array of periods - /// Array of EMA series - public static TSeries[] Calculate(TSeries source, int[] periods) - { - var emaVector = new EmaVector(periods); - return emaVector.Calculate(source); - } -} diff --git a/lib/averages/sma/Sma.Notebook.dib b/lib/averages/sma/Sma.Notebook.dib index fd4ff04d..1b0529a2 100644 --- a/lib/averages/sma/Sma.Notebook.dib +++ b/lib/averages/sma/Sma.Notebook.dib @@ -21,7 +21,8 @@ This notebook demonstrates: 1. **Manual Data Processing**: Understanding Batch vs. Streaming modes. 2. **Streaming with `isNew`**: Handling intra-bar updates. 3. **Large Dataset Processing**: Using Geometric Brownian Motion (GBM) generated data. -4. **Vectorized Operations**: Calculating multiple SMAs simultaneously. +4. **Handling Invalid Values**: Last-value substitution for NaN/Infinity. +5. **SMA vs EMA**: Comparing Simple and Exponential Moving Averages. #!csharp @@ -189,63 +190,9 @@ Console.WriteLine($"Match: {Math.Abs(batchLargeResult.Last().Value - lastStreamV #!markdown -## 4. Vectorized SMA (Multiple Periods) +## 4. Handling Invalid Values (NaN/Infinity) -`SmaVector` allows calculating multiple SMAs (e.g., 5, 10, 20) simultaneously. This is useful for comparing different timeframes. - -### Vectorized Batch - -#!csharp - -int[] periods = { 5, 10, 20 }; -Console.WriteLine($"\n--- Vectorized Batch SMA (Periods: {string.Join(", ", periods)}) ---"); - -var smaVectorBatch = new SmaVector(periods); -var vectorBatchResults = smaVectorBatch.Calculate(closeSeries); - -for (int i = 0; i < periods.Length; i++) -{ - Console.WriteLine($"SMA({periods[i]}) Last Value: {vectorBatchResults[i].Last().Value:F2}"); -} - -#!markdown - -### Vectorized Streaming - -#!csharp - -Console.WriteLine($"\n--- Vectorized Streaming SMA (Periods: {string.Join(", ", periods)}) ---"); - -var smaVectorStream = new SmaVector(periods); -TValue[] lastVectorVal = null; - -foreach(var item in closeSeries) -{ - lastVectorVal = smaVectorStream.Update(item); -} - -for (int i = 0; i < periods.Length; i++) -{ - Console.WriteLine($"SMA({periods[i]}) Last Value: {lastVectorVal[i].Value:F2}"); -} - -// Verification -bool allMatch = true; -for (int i = 0; i < periods.Length; i++) -{ - if (Math.Abs(vectorBatchResults[i].Last().Value - lastVectorVal[i].Value) > 1e-10) - { - allMatch = false; - break; - } -} -Console.WriteLine($"\nAll Vectorized Stream/Batch values match: {allMatch}"); - -#!markdown - -## 5. Handling Invalid Values (NaN/Infinity) - -Both `Sma` and `SmaVector` use **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. +`Sma` uses **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. #!csharp @@ -296,30 +243,9 @@ for (int i = 0; i < seriesWithNaN.Count; i++) Console.WriteLine($" {inputStr,-10} → {output:F2} (IsFinite: {double.IsFinite(output)})"); } -#!csharp - -Console.WriteLine("\n--- Vectorized SMA with Invalid Values ---"); - -int[] periodsNaN = { 5, 10 }; -var smaVectorNaN = new SmaVector(periodsNaN); - -// Feed values including invalid ones -var inputsNaN = new double[] { 100, 110, double.NaN, 120, double.PositiveInfinity, 130 }; -var time = DateTime.Now; - -foreach (var val in inputsNaN) -{ - var results = smaVectorNaN.Update(new TValue(time, val)); - var inputStr = double.IsFinite(val) ? val.ToString("F2") : val.ToString(); - Console.WriteLine($"Input: {inputStr,-10} → SMA(5): {results[0].Value:F2}, SMA(10): {results[1].Value:F2}"); - time = time.AddMinutes(1); -} - -Console.WriteLine("\nAll outputs are finite - invalid inputs were substituted with last valid values."); - #!markdown -## 6. SMA vs EMA Comparison +## 5. SMA vs EMA Comparison The SMA and EMA are both trend-following indicators, but they weight data differently: diff --git a/lib/averages/sma/Sma.cs b/lib/averages/sma/Sma.cs index 0a585869..4b394357 100644 --- a/lib/averages/sma/Sma.cs +++ b/lib/averages/sma/Sma.cs @@ -38,15 +38,13 @@ public sealed class Sma private readonly int _period; private readonly RingBuffer _buffer; - // Running sum maintained separately for O(1) bar correction private double _sum; - private double _p_sum; // Sum AFTER last isNew=true (for correction restore) - private double _p_lastInput; // Input that was added on last isNew=true + private double _p_sum; + private double _p_lastInput; private double _lastValidValue; private double _p_lastValidValue; - private int _tickCount; // Counter for periodic sum resync + private int _tickCount; - // Resync interval: recalculate sum from buffer every N ticks to prevent drift private const int ResyncInterval = 1000; /// @@ -93,23 +91,15 @@ public sealed class Sma return _lastValidValue; } - /// - /// Updates internal state with a new value. - /// Shared logic for both streaming and batch-reconstruction. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private void UpdateState(double val) { - // Calculate what to remove from sum (oldest value if buffer full) double removedValue = _buffer.Count == _buffer.Capacity ? _buffer.Oldest : 0.0; - // Update sum: remove oldest, add newest _sum = _sum - removedValue + val; - // Update buffer _buffer.Add(val); - // Periodic resync: recalculate sum from scratch to eliminate floating-point drift _tickCount++; if (_buffer.IsFull && _tickCount >= ResyncInterval) { @@ -118,43 +108,27 @@ public sealed class Sma } } - /// - /// Updates SMA with the given value. - /// O(1) for both isNew=true and isNew=false. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Current SMA value [MethodImpl(MethodImplOptions.AggressiveInlining)] public TValue Update(TValue input, bool isNew = true) { if (isNew) { - // Get valid value (this may update _lastValidValue) double val = GetValidValue(input.Value); UpdateState(val); - // Save state AFTER this update for potential future corrections _p_sum = _sum; _p_lastInput = val; _p_lastValidValue = _lastValidValue; } else { - // Bar correction: restore to state AFTER last isNew=true, then swap last value - // Restore _lastValidValue BEFORE calling GetValidValue _lastValidValue = _p_lastValidValue; - // Get valid value (this may update _lastValidValue) double val = GetValidValue(input.Value); - // _p_sum is the sum AFTER the last isNew=true completed - // _p_lastInput is the value that was added on last isNew=true - // We want: new_sum = _p_sum - _p_lastInput + val _sum = _p_sum - _p_lastInput + val; - // Update buffer's newest value _buffer.UpdateNewest(val); } @@ -163,11 +137,6 @@ public sealed class Sma return Value; } - /// - /// Updates SMA with the entire series. - /// - /// Input series - /// SMA series public TSeries Update(TSeries source) { if (source.Count == 0) return new TSeries(new List(), new List()); @@ -183,25 +152,15 @@ public sealed class Sma var sourceValues = source.Values; var sourceTimes = source.Times; - // 1. Fast Batch Calculation (SIMD optimized) Calculate(sourceValues, vSpan, _period); - // 2. Copy Times sourceTimes.CopyTo(tSpan); - // 3. Reconstruct State for subsequent updates - // We need to restore _buffer, _sum, and _lastValidValue to what they would be - // if we had processed the series sequentially. - - // Find the last valid value before the reconstruction window - // The reconstruction window is the last 'period' elements (or less if len < period) int windowSize = Math.Min(len, _period); int startIndex = len - windowSize; - // Restore _lastValidValue from before the window if (startIndex > 0) { - // Scan backwards to find last valid value for (int i = startIndex - 1; i >= 0; i--) { if (double.IsFinite(sourceValues[i])) @@ -213,10 +172,9 @@ public sealed class Sma } else { - _lastValidValue = 0; // Reset if starting from 0 + _lastValidValue = 0; } - // Rebuild buffer and sum from last 'period' values using shared logic _buffer.Clear(); _sum = 0; _tickCount = 0; @@ -227,7 +185,6 @@ public sealed class Sma UpdateState(val); } - // Save state for potential future corrections _p_sum = _sum; _p_lastInput = sourceValues[len - 1]; _p_lastValidValue = _lastValidValue; @@ -281,17 +238,11 @@ public sealed class Sma CalculateScalarCore(source, output, period); } - /// - /// Scalar implementation with NaN handling via last-value substitution. - /// Uses circular buffer for sliding window calculation. - /// Optimized with split loops and periodic resync. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private static void CalculateScalarCore(ReadOnlySpan source, Span output, int period) { int len = source.Length; - // Use stackalloc for small periods, otherwise fall back to heap allocation const int StackAllocThreshold = 256; Span buffer = period <= StackAllocThreshold ? stackalloc double[period] @@ -302,8 +253,6 @@ public sealed class Sma int bufferIndex = 0; int i = 0; - // Phase 1: Warmup (0 to period-1) - // No need to remove oldest value, just accumulate int warmupEnd = Math.Min(period, len); for (; i < warmupEnd; i++) { @@ -318,9 +267,6 @@ public sealed class Sma output[i] = sum / (i + 1); } - // Phase 2: Hot loop (period to len) - // Buffer is full, remove oldest, add newest - // Optimized buffer indexing (no modulo) int tickCount = 0; for (; i < len; i++) { @@ -330,23 +276,19 @@ public sealed class Sma else val = lastValid; - // Remove oldest, add newest sum = sum - buffer[bufferIndex] + val; buffer[bufferIndex] = val; - // Increment buffer index with wrap-around check (faster than modulo) bufferIndex++; if (bufferIndex >= period) bufferIndex = 0; output[i] = sum / period; - // Periodic resync every 1000 ticks tickCount++; if (tickCount >= ResyncInterval) { tickCount = 0; - // Recalculate sum from buffer to prevent drift double recalcSum = 0; for (int k = 0; k < period; k++) { @@ -357,29 +299,17 @@ public sealed class Sma } } - /// - /// SIMD-optimized implementation for SMA calculation. - /// Processes 4 consecutive values per iteration using AVX2 (Vector256<double>). - /// Assumes input contains no NaN/Infinity values. - /// - /// - /// Key insight: For consecutive positions i, i+1, i+2, i+3: - /// - sum[i+1] = sum[i] - src[i-period+1] + src[i+1] - /// - We can vectorize the load of 4 "leaving" values and 4 "entering" values - /// - Then use prefix-sum style to compute the 4 sums from one base sum - /// [MethodImpl(MethodImplOptions.AggressiveOptimization)] private static unsafe void CalculateSimdCore(ReadOnlySpan source, Span output, int period) { int len = source.Length; - const int VectorWidth = 4; // Vector256 holds 4 doubles + const int VectorWidth = 4; fixed (double* srcPtr = source) fixed (double* outPtr = output) { double invPeriod = 1.0 / period; - // Phase 1: Warmup - scalar processing until buffer is full int warmupEnd = Math.Min(period, len); double sum = 0; for (int i = 0; i < warmupEnd; i++) @@ -391,8 +321,6 @@ public sealed class Sma if (len <= period) return; - // Phase 2: SIMD hot loop - // Uses prefix-sum approach to break dependency chain var vInvPeriod = Vector256.Create(invPeriod); var vZero = Vector256.Zero; int simdEnd = period + ((len - period) / VectorWidth) * VectorWidth; @@ -400,44 +328,31 @@ public sealed class Sma for (int i = period; i < simdEnd; i += VectorWidth) { - // Load 4 entering values and 4 leaving values var vNew = Avx.LoadVector256(srcPtr + i); var vOld = Avx.LoadVector256(srcPtr + i - period); - // Delta = New - Old var vDelta = Avx.Subtract(vNew, vOld); - // Prefix sum of Deltas - // Step 1: Shift right by 1 element (insert 0) - // [D0, D1, D2, D3] -> [0, D0, D1, D2] var vShift1 = Avx2.Permute4x64(vDelta.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShift1 = Avx.Blend(vZero, vShift1, 0b_1110); - var vP1 = Avx.Add(vDelta, vShift1); // [D0, D0+D1, D1+D2, D2+D3] + var vP1 = Avx.Add(vDelta, vShift1); - // Step 2: Shift right by 2 elements (insert 0) - // [D0, D0+D1, D1+D2, D2+D3] -> [0, 0, D0, D0+D1] var vShift2 = Avx2.Permute4x64(vP1.AsUInt64(), 0b_01_00_00_00).AsDouble(); vShift2 = Avx.Blend(vZero, vShift2, 0b_1100); - var vP2 = Avx.Add(vP1, vShift2); // [D0, D0+D1, D0+D1+D2, D0+D1+D2+D3] + var vP2 = Avx.Add(vP1, vShift2); - // Add previous sum to all var vSumPrev = Vector256.Create(sum); var vSums = Avx.Add(vSumPrev, vP2); - // Store result var vResult = Avx.Multiply(vSums, vInvPeriod); Avx.Store(outPtr + i, vResult); - // Update sum for next iteration (last element of vSums) sum = vSums.GetElement(3); - // Periodic resync every 1000 ticks tickCount += VectorWidth; if (tickCount >= ResyncInterval) { tickCount = 0; - // Recalculate sum from scratch using the window ending at i + VectorWidth - 1 - // Window: [i + VectorWidth - period ... i + VectorWidth - 1] int lastIdx = i + VectorWidth - 1; double recalcSum = 0; for (int k = 0; k < period; k++) @@ -448,7 +363,6 @@ public sealed class Sma } } - // Phase 3: Scalar tail for (int i = simdEnd; i < len; i++) { sum = sum - srcPtr[i - period] + srcPtr[i]; diff --git a/lib/averages/sma/Sma.md b/lib/averages/sma/Sma.md index c000fee4..5991fe6a 100644 --- a/lib/averages/sma/Sma.md +++ b/lib/averages/sma/Sma.md @@ -99,33 +99,9 @@ Console.WriteLine($"Last SMA: {smaOutput[^1]}"); * **2-3x faster** than TSeries API for large datasets * **Compatible** with `ArrayPool` for buffer management -### Multi-Period SMA (`SmaVector`) - -The `SmaVector` class calculates multiple SMAs with different periods on the same input series simultaneously. - -```csharp -using QuanTAlib; - -// Initialize with multiple periods -int[] periods = { 5, 10, 20 }; -var smaVector = new SmaVector(periods); - -// Streaming update -TValue[] results = smaVector.Update(new TValue(time, price)); - -// Access values -Console.WriteLine($"SMA(5): {results[0].Value}"); -Console.WriteLine($"SMA(10): {results[1].Value}"); -Console.WriteLine($"SMA(20): {results[2].Value}"); - -// Batch calculation -TSeries source = ...; -TSeries[] seriesResults = smaVector.Calculate(source); -``` - ### Bar Correction (isNew Parameter) -Both `Sma` and `SmaVector` support intra-bar updates for real-time trading systems: +`Sma` supports intra-bar updates for real-time trading systems: ```csharp var sma = new Sma(10); @@ -151,7 +127,7 @@ sma.Update(new TValue(time + 1, 101.2), isNew: true); ### Handling Invalid Values (NaN/Infinity) -Both `Sma` and `SmaVector` use **last-value substitution** for handling invalid inputs: +`Sma` uses **last-value substitution** for handling invalid inputs: ```csharp var sma = new Sma(10); diff --git a/lib/averages/sma/SmaVector.Tests.cs b/lib/averages/sma/SmaVector.Tests.cs deleted file mode 100644 index 99209d59..00000000 --- a/lib/averages/sma/SmaVector.Tests.cs +++ /dev/null @@ -1,371 +0,0 @@ - -namespace QuanTAlib.Tests; - -public class SmaVectorTests -{ - [Fact] - public void Initialization_WithPeriods_Works() - { - int[] periods = { 5, 10, 20 }; - var smaVector = new SmaVector(periods); - - var res = smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(3, res.Length); - Assert.Equal(100.0, res[0].Value, 1e-9); - Assert.Equal(100.0, res[1].Value, 1e-9); - Assert.Equal(100.0, res[2].Value, 1e-9); - } - - [Fact] - public void Initialization_WithZeroPeriod_ThrowsArgumentException() - { - int[] periods = { 10, 0, 20 }; - - Assert.Throws(() => new SmaVector(periods)); - } - - [Fact] - public void Initialization_WithNegativePeriod_ThrowsArgumentException() - { - int[] periods = { 10, -5, 20 }; - - Assert.Throws(() => new SmaVector(periods)); - } - - [Fact] - public void Calc_Streaming_MatchesSingleSma() - { - int[] periods = { 5, 10, 20 }; - var smaVector = new SmaVector(periods); - var smaSingles = periods.Select(p => new Sma(p)).ToArray(); - - var values = new double[] { 10, 20, 30, 40, 50, 40, 30, 20, 10 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = smaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = smaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - Assert.Equal(singleRes.Time, multiRes[i].Time); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Calc_Series_MatchesSingleSma() - { - int[] periods = { 5, 10, 20 }; - var smaVector = new SmaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var multiRes = smaVector.Calculate(series); - - // Reset and recalculate for comparison - var smaSingles = periods.Select(p => new Sma(p)).ToArray(); - for (int j = 0; j < len; j++) - { - var tVal = new TValue(new DateTime(t[j], DateTimeKind.Utc), v[j]); - for (int i = 0; i < periods.Length; i++) - { - var singleRes = smaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Values[j], 1e-8); - } - } - } - - [Fact] - public void Calc_Series_MatchesStreaming() - { - int[] periods = { 5, 10, 20 }; - var smaVectorBatch = new SmaVector(periods); - var smaVectorStream = new SmaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var batchRes = smaVectorBatch.Calculate(series); - - for (int i = 0; i < len; i++) - { - var tVal = new TValue(new DateTime(t[i], DateTimeKind.Utc), v[i]); - var streamRes = smaVectorStream.Update(tVal); - - for (int j = 0; j < periods.Length; j++) - { - Assert.Equal(batchRes[j].Values[i], streamRes[j].Value, 1e-9); - } - } - } - - [Fact] - public void Calculate_Static_MatchesInstanceMethod() - { - int[] periods = { 5, 10, 20 }; - - int len = 50; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var instanceSma = new SmaVector(periods); - var instanceRes = instanceSma.Calculate(series); - - var staticRes = SmaVector.Calculate(series, periods); - - for (int i = 0; i < periods.Length; i++) - { - Assert.Equal(instanceRes[i].Count, staticRes[i].Count); - for (int j = 0; j < len; j++) - { - Assert.Equal(instanceRes[i].Values[j], staticRes[i].Values[j], 1e-9); - } - } - } - - [Fact] - public void Reset_ClearsState() - { - int[] periods = { 10 }; - var smaVector = new SmaVector(periods); - - smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - smaVector.Update(new TValue(DateTime.UtcNow, 200.0)); - smaVector.Reset(); - - var res = smaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - - Assert.Equal(50.0, res[0].Value, 1e-9); - } - - [Fact] - public void Update_NaN_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var smaVector = new SmaVector(periods); - - smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - smaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterNaN = smaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in resultAfterNaN) - { - Assert.True(double.IsFinite(result.Value), $"Expected finite value but got {result.Value}"); - } - } - - [Fact] - public void Update_Infinity_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var smaVector = new SmaVector(periods); - - smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - smaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterPosInf = smaVector.Update(new TValue(DateTime.UtcNow, double.PositiveInfinity)); - foreach (var result in resultAfterPosInf) - { - Assert.True(double.IsFinite(result.Value)); - } - - var resultAfterNegInf = smaVector.Update(new TValue(DateTime.UtcNow, double.NegativeInfinity)); - foreach (var result in resultAfterNegInf) - { - Assert.True(double.IsFinite(result.Value)); - } - } - - [Fact] - public void Update_MultipleNaN_ContinuesWithLastValid() - { - int[] periods = { 5, 10 }; - var smaVector = new SmaVector(periods); - - smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - smaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - smaVector.Update(new TValue(DateTime.UtcNow, 120.0)); - - var r1 = smaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r2 = smaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r3 = smaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in r1) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r2) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r3) Assert.True(double.IsFinite(result.Value)); - } - - [Fact] - public void Calculate_Series_HandlesNaN() - { - int[] periods = { 5, 10 }; - var smaVector = new SmaVector(periods); - - var t = new System.Collections.Generic.List(); - var v = new System.Collections.Generic.List(); - var now = DateTime.UtcNow; - - t.Add(now.Ticks); v.Add(100.0); - t.Add(now.AddMinutes(1).Ticks); v.Add(110.0); - t.Add(now.AddMinutes(2).Ticks); v.Add(double.NaN); - t.Add(now.AddMinutes(3).Ticks); v.Add(120.0); - t.Add(now.AddMinutes(4).Ticks); v.Add(double.PositiveInfinity); - t.Add(now.AddMinutes(5).Ticks); v.Add(130.0); - - var series = new TSeries(t, v); - var results = smaVector.Calculate(series); - - foreach (var periodResults in results) - { - foreach (var val in periodResults.Values) - { - Assert.True(double.IsFinite(val), $"Expected finite value but got {val}"); - } - } - } - - [Fact] - public void Reset_ClearsLastValidValue() - { - int[] periods = { 10 }; - var smaVector = new SmaVector(periods); - - smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - smaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - smaVector.Reset(); - - var result = smaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - Assert.Equal(50.0, result[0].Value, 1e-9); - } - - [Fact] - public void NaN_Handling_MatchesSingleSma() - { - int[] periods = { 5, 10, 20 }; - var smaVector = new SmaVector(periods); - var smaSingles = periods.Select(p => new Sma(p)).ToArray(); - - var values = new double[] { 10, 20, double.NaN, 40, double.PositiveInfinity, 60, 70 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = smaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = smaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Values_Property_UpdatesAfterUpdate() - { - int[] periods = { 5, 10 }; - var smaVector = new SmaVector(periods); - - var result = smaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(result[0].Value, smaVector.Values[0].Value); - Assert.Equal(result[1].Value, smaVector.Values[1].Value); - } - - [Fact] - public void Values_Property_UpdatesAfterCalculate() - { - int[] periods = { 5, 10 }; - var smaVector = new SmaVector(periods); - - var t = new System.Collections.Generic.List { 100, 200, 300 }; - var v = new System.Collections.Generic.List { 10.0, 20.0, 30.0 }; - var series = new TSeries(t, v); - - var results = smaVector.Calculate(series); - - Assert.Equal(results[0].Last.Value, smaVector.Values[0].Value, 1e-9); - Assert.Equal(results[1].Last.Value, smaVector.Values[1].Value, 1e-9); - } - - [Fact] - public void Update_BarCorrection_WorksCorrectly() - { - int[] periods = { 3 }; - var smaVector = new SmaVector(periods); - - smaVector.Update(new TValue(DateTime.UtcNow, 10.0), isNew: true); - smaVector.Update(new TValue(DateTime.UtcNow, 20.0), isNew: true); - smaVector.Update(new TValue(DateTime.UtcNow, 30.0), isNew: true); - - var res1 = smaVector.Values[0].Value; - Assert.Equal(20.0, res1, 1e-9); // (10+20+30)/3 = 20 - - // Correct the last bar - var res2 = smaVector.Update(new TValue(DateTime.UtcNow, 60.0), isNew: false); - - Assert.Equal(30.0, res2[0].Value, 1e-9); // (10+20+60)/3 = 30 - } - - [Fact] - public void SMA_MatchesExpectedValues() - { - int[] periods = { 3 }; - var smaVector = new SmaVector(periods); - - // Test sequence: 10, 20, 30, 40, 50 - // Expected SMA(3): 10, 15, 20, 30, 40 - var expected = new double[] { 10, 15, 20, 30, 40 }; - var values = new double[] { 10, 20, 30, 40, 50 }; - var time = DateTime.UtcNow; - - for (int i = 0; i < values.Length; i++) - { - var res = smaVector.Update(new TValue(time, values[i])); - Assert.Equal(expected[i], res[0].Value, 1e-9); - time = time.AddMinutes(1); - } - } -} diff --git a/lib/averages/sma/SmaVector.cs b/lib/averages/sma/SmaVector.cs deleted file mode 100644 index 7e0ac6da..00000000 --- a/lib/averages/sma/SmaVector.cs +++ /dev/null @@ -1,185 +0,0 @@ -using System.Runtime.CompilerServices; -using System.Runtime.InteropServices; - -namespace QuanTAlib; - -/// -/// Multi-Period Simple Moving Average (SMA) - SIMD optimized. -/// Calculates multiple SMAs with different periods for the same input series in parallel. -/// Uses last-value substitution for invalid inputs (NaN/Infinity). -/// -[SkipLocalsInit] -public class SmaVector -{ - private readonly RingBuffer[] _buffers; - private readonly RingBuffer[] _p_buffers; // Previous state for bar correction - private readonly int _count; - private double _lastValidValue; - - /// - /// Current SMA values for all periods. - /// - public ReadOnlySpan Values => _values; - - private readonly TValue[] _values; - - /// - /// Initializes SmaVector with specified periods. - /// - /// Array of periods (each must be > 0) - public SmaVector(int[] periods) - { - _count = periods.Length; - _buffers = new RingBuffer[_count]; - _p_buffers = new RingBuffer[_count]; - _values = new TValue[_count]; - - for (int i = 0; i < _count; i++) - { - ArgumentOutOfRangeException.ThrowIfLessThanOrEqual(periods[i], 0); - _buffers[i] = new RingBuffer(periods[i]); - _p_buffers[i] = new RingBuffer(periods[i]); - } - } - - /// - /// Gets a valid input value, using last-value substitution for non-finite inputs. - /// - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private double GetValidValue(double input) - { - if (double.IsFinite(input)) - { - _lastValidValue = input; - return input; - } - return _lastValidValue; - } - - /// - /// Resets all SMA states. - /// - public void Reset() - { - for (int i = 0; i < _count; i++) - { - _buffers[i].Clear(); - _p_buffers[i].Clear(); - } - _lastValidValue = 0; - Array.Clear(_values); - } - - /// - /// Updates SMAs with the given value. - /// Uses last-value substitution: invalid inputs (NaN/Infinity) are replaced with - /// the last known good value, providing continuity in the output series. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Array of SMA values - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public TValue[] Update(TValue input, bool isNew = true) - { - if (isNew) - { - // Save current state for potential bar correction - for (int i = 0; i < _count; i++) - { - _p_buffers[i].CopyFrom(_buffers[i]); - } - } - else - { - // Restore previous state for bar correction - for (int i = 0; i < _count; i++) - { - _buffers[i].CopyFrom(_p_buffers[i]); - } - } - - // Last-value substitution: replace non-finite inputs with last valid value - double val = GetValidValue(input.Value); - - // Update each buffer and calculate SMA - for (int i = 0; i < _count; i++) - { - _buffers[i].Add(val); - _values[i] = new TValue(input.Time, _buffers[i].Average); - } - - return _values; - } - - /// - /// Calculates SMAs for the entire series. - /// - /// Input series - /// Array of SMA series - public TSeries[] Calculate(TSeries source) - { - int len = source.Count; - var resultSeries = new TSeries[_count]; - - // Reset state for fresh calculation - for (int i = 0; i < _count; i++) - { - _buffers[i].Clear(); - } - _lastValidValue = 0; - - // Pre-allocate lists - var tLists = new List[_count]; - var vLists = new List[_count]; - - for (int i = 0; i < _count; i++) - { - tLists[i] = new List(len); - vLists[i] = new List(len); - CollectionsMarshal.SetCount(tLists[i], len); - CollectionsMarshal.SetCount(vLists[i], len); - } - - var sourceValues = source.Values; - var sourceTimes = source.Times; - - for (int t = 0; t < len; t++) - { - double val = sourceValues[t]; - long time = sourceTimes[t]; - - // Last-value substitution: replace non-finite inputs with last valid value - val = GetValidValue(val); - - for (int i = 0; i < _count; i++) - { - _buffers[i].Add(val); - CollectionsMarshal.AsSpan(tLists[i])[t] = time; - CollectionsMarshal.AsSpan(vLists[i])[t] = _buffers[i].Average; - } - } - - // Create TSeries and update Values - for (int i = 0; i < _count; i++) - { - resultSeries[i] = new TSeries(tLists[i], vLists[i]); - var lastT = CollectionsMarshal.AsSpan(tLists[i])[len - 1]; - var lastV = CollectionsMarshal.AsSpan(vLists[i])[len - 1]; - _values[i] = new TValue(lastT, lastV); - } - - return resultSeries; - } - - /// - /// Calculates SMAs for the entire series using specified periods. - /// - /// Input series - /// Array of periods - /// Array of SMA series - public static TSeries[] Calculate(TSeries source, int[] periods) - { - var smaVector = new SmaVector(periods); - return smaVector.Calculate(source); - } -} diff --git a/lib/averages/trima/Trima.Tests.cs b/lib/averages/trima/Trima.Tests.cs index 1d7158ed..82cbcaa2 100644 --- a/lib/averages/trima/Trima.Tests.cs +++ b/lib/averages/trima/Trima.Tests.cs @@ -108,7 +108,6 @@ public class TrimaTests trima.Update(new TValue(DateTime.UtcNow, 100)); trima.Update(new TValue(DateTime.UtcNow, 105)); - double valueBefore = trima.Value; trima.Reset(); @@ -151,20 +150,24 @@ public class TrimaTests // Calculate iteratively var iterativeResults = new TSeries(); +#pragma warning disable S4158 // Collection is known to be empty foreach (var item in series) { iterativeResults.Add(trimaIterative.Update(item)); } +#pragma warning restore S4158 // Calculate batch var batchResults = trimaBatch.Update(series); // Compare Assert.Equal(iterativeResults.Count, batchResults.Count); +#pragma warning disable S2583 // Condition always evaluates to false for (int i = 0; i < iterativeResults.Count; i++) { Assert.Equal(iterativeResults[i].Value, batchResults[i].Value, 1e-10); } +#pragma warning restore S2583 } [Fact] diff --git a/lib/averages/trima/Trima.cs b/lib/averages/trima/Trima.cs index 5c29d4e3..a72b288d 100644 --- a/lib/averages/trima/Trima.cs +++ b/lib/averages/trima/Trima.cs @@ -8,29 +8,17 @@ namespace QuanTAlib; /// TRIMA: Triangular Moving Average /// /// -/// TRIMA is a weighted moving average where the weights increase linearly to the middle -/// of the period and then decrease linearly. It places the most weight on the middle -/// portion of the data series. +/// TRIMA is a weighted moving average where weights increase linearly to the middle +/// and then decrease. It is equivalent to a double SMA: SMA(SMA(period1), period2). /// /// Calculation: -/// TRIMA(period) = SMA(SMA(period1), period2) -/// where: /// period1 = period / 2 + 1 /// period2 = (period + 1) / 2 /// -/// This implementation uses a flattened structure with two internal SMA buffers -/// to ensure correct handling of warmup periods and bar corrections without -/// the overhead of composed objects. -/// -/// Key characteristics: -/// - Smoother than SMA -/// - Double smoothing (lag is higher than SMA) -/// - Weights form a triangle +/// Characteristics: +/// - Smoother than SMA, higher lag /// - O(1) time complexity /// - O(period) space complexity -/// -/// Sources: -/// - https://www.investopedia.com/terms/t/triangularaverage.asp /// [SkipLocalsInit] public sealed class Trima @@ -41,36 +29,22 @@ public sealed class Trima private readonly RingBuffer _buffer1; private readonly RingBuffer _buffer2; - // SMA1 State - private double _sum1; - private double _p_sum1; - private double _p_lastInput1; - private double _lastValidValue1; - private double _p_lastValidValue1; + private double _sum1, _p_sum1, _p_lastInput1, _lastValidValue1, _p_lastValidValue1; private int _tickCount1; - // SMA2 State - private double _sum2; - private double _p_sum2; - private double _p_lastInput2; + private double _sum2, _p_sum2, _p_lastInput2; private int _tickCount2; private int _sampleCount; private const int ResyncInterval = 1000; - /// - /// Display name for the indicator. - /// public string Name { get; } + public TValue Value { get; private set; } + public bool IsHot => _sampleCount >= _period; - /// - /// Creates TRIMA with specified period. - /// - /// Number of values to average (must be > 0) public Trima(int period) { - if (period <= 0) - throw new ArgumentException("Period must be greater than 0", nameof(period)); + if (period <= 0) throw new ArgumentException("Period must be greater than 0", nameof(period)); _period = period; _p1 = period / 2 + 1; @@ -82,19 +56,6 @@ public sealed class Trima Name = $"Trima({period})"; } - /// - /// Current TRIMA value. - /// - public TValue Value { get; private set; } - - /// - /// True if the TRIMA has enough data to produce valid results. - /// - public bool IsHot => _sampleCount >= _period; - - /// - /// Gets a valid input value, using last-value substitution for non-finite inputs. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private double GetValidValue(double input) { @@ -106,12 +67,6 @@ public sealed class Trima return _lastValidValue1; } - /// - /// Updates TRIMA with the given value. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Current TRIMA value [MethodImpl(MethodImplOptions.AggressiveInlining)] public TValue Update(TValue input, bool isNew = true) { @@ -119,13 +74,12 @@ public sealed class Trima { _sampleCount++; - // SMA 1 Update + // SMA 1 double val1 = GetValidValue(input.Value); double removed1 = _buffer1.Count == _buffer1.Capacity ? _buffer1.Oldest : 0.0; _sum1 = _sum1 - removed1 + val1; _buffer1.Add(val1); - // Resync SMA1 _tickCount1++; if (_buffer1.IsFull && _tickCount1 >= ResyncInterval) { @@ -133,21 +87,17 @@ public sealed class Trima _sum1 = _buffer1.Sum(); } - // Save SMA1 state _p_sum1 = _sum1; _p_lastInput1 = val1; _p_lastValidValue1 = _lastValidValue1; - // SMA 1 Result double sma1Result = _sum1 / _buffer1.Count; - // SMA 2 Update (Input is sma1Result) - // Note: sma1Result is always finite if input stream has at least one finite value + // SMA 2 double removed2 = _buffer2.Count == _buffer2.Capacity ? _buffer2.Oldest : 0.0; _sum2 = _sum2 - removed2 + sma1Result; _buffer2.Add(sma1Result); - // Resync SMA2 _tickCount2++; if (_buffer2.IsFull && _tickCount2 >= ResyncInterval) { @@ -155,13 +105,10 @@ public sealed class Trima _sum2 = _buffer2.Sum(); } - // Save SMA2 state _p_sum2 = _sum2; _p_lastInput2 = sma1Result; - // Final Result - double trimaResult = _sum2 / _buffer2.Count; - Value = new TValue(input.Time, trimaResult); + Value = new TValue(input.Time, _sum2 / _buffer2.Count); } else { @@ -177,23 +124,16 @@ public sealed class Trima _sum2 = _p_sum2 - _p_lastInput2 + sma1Result; _buffer2.UpdateNewest(sma1Result); - double trimaResult = _sum2 / _buffer2.Count; - Value = new TValue(input.Time, trimaResult); + Value = new TValue(input.Time, _sum2 / _buffer2.Count); } return Value; } - /// - /// Updates TRIMA with the entire series. - /// - /// Input series - /// TRIMA series public TSeries Update(TSeries source) { if (source.Count == 0) return new TSeries(new List(), new List()); - // Use the static Calculate method for performance int len = source.Count; var t = new List(len); var v = new List(len); @@ -202,43 +142,30 @@ public sealed class Trima var tSpan = CollectionsMarshal.AsSpan(t); var vSpan = CollectionsMarshal.AsSpan(v); - var sourceValues = source.Values; - var sourceTimes = source.Times; + + Calculate(source.Values, vSpan, _period); + source.Times.CopyTo(tSpan); - Calculate(sourceValues, vSpan, _period); - sourceTimes.CopyTo(tSpan); - - // Restore state by replaying the last part - // We need to replay enough to fill both SMAs + // Restore state int lookback = _p1 + _p2; int startIndex = Math.Max(0, len - lookback); - - // Reset internal state Reset(); - // Replay for (int i = startIndex; i < len; i++) { - Update(new TValue(sourceTimes[i], sourceValues[i]), isNew: true); + Update(new TValue(source.Times[i], source.Values[i]), isNew: true); } Value = new TValue(tSpan[len - 1], vSpan[len - 1]); return new TSeries(t, v); } - /// - /// Calculates TRIMA for the entire series using a new instance. - /// public static TSeries Calculate(TSeries source, int period) { var trima = new Trima(period); return trima.Update(source); } - /// - /// Calculates TRIMA in-place. - /// Uses ArrayPool to allocate temporary buffer and chains optimized SMA calculations. - /// public static void Calculate(ReadOnlySpan source, Span output, int period) { if (source.Length != output.Length) @@ -249,16 +176,12 @@ public sealed class Trima int p1 = period / 2 + 1; int p2 = (period + 1) / 2; - // Rent a temporary buffer for the intermediate SMA double[] tempArray = ArrayPool.Shared.Rent(source.Length); Span tempSpan = tempArray.AsSpan(0, source.Length); try { - // SMA 1 Sma.Calculate(source, tempSpan, p1); - - // SMA 2 (TRIMA) Sma.Calculate(tempSpan, output, p2); } finally @@ -267,24 +190,15 @@ public sealed class Trima } } - /// - /// Resets the TRIMA state. - /// public void Reset() { _buffer1.Clear(); _buffer2.Clear(); - _sum1 = 0; - _p_sum1 = 0; - _p_lastInput1 = 0; - _lastValidValue1 = 0; - _p_lastValidValue1 = 0; + _sum1 = _p_sum1 = _p_lastInput1 = _lastValidValue1 = _p_lastValidValue1 = 0; _tickCount1 = 0; - _sum2 = 0; - _p_sum2 = 0; - _p_lastInput2 = 0; + _sum2 = _p_sum2 = _p_lastInput2 = 0; _tickCount2 = 0; _sampleCount = 0; diff --git a/lib/averages/trima/TrimaVector.Tests.cs b/lib/averages/trima/TrimaVector.Tests.cs deleted file mode 100644 index 63d1ff66..00000000 --- a/lib/averages/trima/TrimaVector.Tests.cs +++ /dev/null @@ -1,363 +0,0 @@ -namespace QuanTAlib.Tests; - -public class TrimaVectorTests -{ - [Fact] - public void Initialization_WithPeriods_Works() - { - int[] periods = { 5, 10, 20 }; - var trimaVector = new TrimaVector(periods); - - var res = trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(3, res.Length); - Assert.Equal(100.0, res[0].Value, 1e-9); - Assert.Equal(100.0, res[1].Value, 1e-9); - Assert.Equal(100.0, res[2].Value, 1e-9); - } - - [Fact] - public void Initialization_WithZeroPeriod_ThrowsArgumentException() - { - int[] periods = { 10, 0, 20 }; - - Assert.Throws(() => new TrimaVector(periods)); - } - - [Fact] - public void Initialization_WithNegativePeriod_ThrowsArgumentException() - { - int[] periods = { 10, -5, 20 }; - - Assert.Throws(() => new TrimaVector(periods)); - } - - [Fact] - public void Calc_Streaming_MatchesSingleTrima() - { - int[] periods = { 5, 10, 20 }; - var trimaVector = new TrimaVector(periods); - var trimaSingles = periods.Select(p => new Trima(p)).ToArray(); - - var values = new double[] { 10, 20, 30, 40, 50, 40, 30, 20, 10 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = trimaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = trimaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - Assert.Equal(singleRes.Time, multiRes[i].Time); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Calc_Series_MatchesSingleTrima() - { - int[] periods = { 5, 10, 20 }; - var trimaVector = new TrimaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var multiRes = trimaVector.Calculate(series); - - // Reset and recalculate for comparison - var trimaSingles = periods.Select(p => new Trima(p)).ToArray(); - for (int j = 0; j < len; j++) - { - var tVal = new TValue(new DateTime(t[j], DateTimeKind.Utc), v[j]); - for (int i = 0; i < periods.Length; i++) - { - var singleRes = trimaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Values[j], 1e-8); - } - } - } - - [Fact] - public void Calc_Series_MatchesStreaming() - { - int[] periods = { 5, 10, 20 }; - var trimaVectorBatch = new TrimaVector(periods); - var trimaVectorStream = new TrimaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var batchRes = trimaVectorBatch.Calculate(series); - - for (int i = 0; i < len; i++) - { - var tVal = new TValue(new DateTime(t[i], DateTimeKind.Utc), v[i]); - var streamRes = trimaVectorStream.Update(tVal); - - for (int j = 0; j < periods.Length; j++) - { - Assert.Equal(batchRes[j].Values[i], streamRes[j].Value, 1e-9); - } - } - } - - [Fact] - public void Calculate_Static_MatchesInstanceMethod() - { - int[] periods = { 5, 10, 20 }; - - int len = 50; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var instanceTrima = new TrimaVector(periods); - var instanceRes = instanceTrima.Calculate(series); - - var staticRes = TrimaVector.Calculate(series, periods); - - for (int i = 0; i < periods.Length; i++) - { - Assert.Equal(instanceRes[i].Count, staticRes[i].Count); - for (int j = 0; j < len; j++) - { - Assert.Equal(instanceRes[i].Values[j], staticRes[i].Values[j], 1e-9); - } - } - } - - [Fact] - public void Reset_ClearsState() - { - int[] periods = { 10 }; - var trimaVector = new TrimaVector(periods); - - trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - trimaVector.Update(new TValue(DateTime.UtcNow, 200.0)); - trimaVector.Reset(); - - var res = trimaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - - Assert.Equal(50.0, res[0].Value, 1e-9); - } - - [Fact] - public void Update_NaN_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var trimaVector = new TrimaVector(periods); - - trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - trimaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterNaN = trimaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in resultAfterNaN) - { - Assert.True(double.IsFinite(result.Value), $"Expected finite value but got {result.Value}"); - } - } - - [Fact] - public void Update_Infinity_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var trimaVector = new TrimaVector(periods); - - trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - trimaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterPosInf = trimaVector.Update(new TValue(DateTime.UtcNow, double.PositiveInfinity)); - foreach (var result in resultAfterPosInf) - { - Assert.True(double.IsFinite(result.Value)); - } - - var resultAfterNegInf = trimaVector.Update(new TValue(DateTime.UtcNow, double.NegativeInfinity)); - foreach (var result in resultAfterNegInf) - { - Assert.True(double.IsFinite(result.Value)); - } - } - - [Fact] - public void Update_MultipleNaN_ContinuesWithLastValid() - { - int[] periods = { 5, 10 }; - var trimaVector = new TrimaVector(periods); - - trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - trimaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - trimaVector.Update(new TValue(DateTime.UtcNow, 120.0)); - - var r1 = trimaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r2 = trimaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r3 = trimaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in r1) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r2) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r3) Assert.True(double.IsFinite(result.Value)); - } - - [Fact] - public void Calculate_Series_HandlesNaN() - { - int[] periods = { 5, 10 }; - var trimaVector = new TrimaVector(periods); - - var t = new System.Collections.Generic.List(); - var v = new System.Collections.Generic.List(); - var now = DateTime.UtcNow; - - t.Add(now.Ticks); v.Add(100.0); - t.Add(now.AddMinutes(1).Ticks); v.Add(110.0); - t.Add(now.AddMinutes(2).Ticks); v.Add(double.NaN); - t.Add(now.AddMinutes(3).Ticks); v.Add(120.0); - t.Add(now.AddMinutes(4).Ticks); v.Add(double.PositiveInfinity); - t.Add(now.AddMinutes(5).Ticks); v.Add(130.0); - - var series = new TSeries(t, v); - var results = trimaVector.Calculate(series); - - foreach (var periodResults in results) - { - foreach (var val in periodResults.Values) - { - Assert.True(double.IsFinite(val), $"Expected finite value but got {val}"); - } - } - } - - [Fact] - public void Reset_ClearsLastValidValue() - { - int[] periods = { 10 }; - var trimaVector = new TrimaVector(periods); - - trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - trimaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - trimaVector.Reset(); - - var result = trimaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - Assert.Equal(50.0, result[0].Value, 1e-9); - } - - [Fact] - public void NaN_Handling_MatchesSingleTrima() - { - int[] periods = { 5, 10, 20 }; - var trimaVector = new TrimaVector(periods); - var trimaSingles = periods.Select(p => new Trima(p)).ToArray(); - - var values = new double[] { 10, 20, double.NaN, 40, double.PositiveInfinity, 60, 70 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = trimaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = trimaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Values_Property_UpdatesAfterUpdate() - { - int[] periods = { 5, 10 }; - var trimaVector = new TrimaVector(periods); - - var result = trimaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(result[0].Value, trimaVector.Values[0].Value); - Assert.Equal(result[1].Value, trimaVector.Values[1].Value); - } - - [Fact] - public void Values_Property_UpdatesAfterCalculate() - { - int[] periods = { 5, 10 }; - var trimaVector = new TrimaVector(periods); - - var t = new System.Collections.Generic.List { 100, 200, 300 }; - var v = new System.Collections.Generic.List { 10.0, 20.0, 30.0 }; - var series = new TSeries(t, v); - - var results = trimaVector.Calculate(series); - - Assert.Equal(results[0].Last.Value, trimaVector.Values[0].Value, 1e-9); - Assert.Equal(results[1].Last.Value, trimaVector.Values[1].Value, 1e-9); - } - - [Fact] - public void Update_BarCorrection_WorksCorrectly() - { - int[] periods = { 3 }; - var trimaVector = new TrimaVector(periods); - - // TRIMA(3) = SMA(SMA(3, 2), 2) - // p1 = 3/2 + 1 = 2 - // p2 = (3+1)/2 = 2 - // SMA1(2): 10 -> 10 - // SMA2(2): 10 -> 10 - trimaVector.Update(new TValue(DateTime.UtcNow, 10.0), isNew: true); - - // SMA1(2): 10, 20 -> 15 - // SMA2(2): 10, 15 -> 12.5 - trimaVector.Update(new TValue(DateTime.UtcNow, 20.0), isNew: true); - - // SMA1(2): 20, 30 -> 25 - // SMA2(2): 15, 25 -> 20 - trimaVector.Update(new TValue(DateTime.UtcNow, 30.0), isNew: true); - - var res1 = trimaVector.Values[0].Value; - Assert.Equal(20.0, res1, 1e-9); - - // Correct the last bar: 30 -> 60 - // SMA1(2): 20, 60 -> 40 - // SMA2(2): 15, 40 -> 27.5 - var res2 = trimaVector.Update(new TValue(DateTime.UtcNow, 60.0), isNew: false); - - Assert.Equal(27.5, res2[0].Value, 1e-9); - } -} diff --git a/lib/averages/trima/TrimaVector.cs b/lib/averages/trima/TrimaVector.cs deleted file mode 100644 index deb4c9b0..00000000 --- a/lib/averages/trima/TrimaVector.cs +++ /dev/null @@ -1,241 +0,0 @@ -using System.Runtime.CompilerServices; -using System.Runtime.InteropServices; - -namespace QuanTAlib; - -/// -/// Multi-Period Triangular Moving Average (TRIMA) - SIMD optimized. -/// Calculates multiple TRIMAs with different periods for the same input series in parallel. -/// Uses last-value substitution for invalid inputs (NaN/Infinity). -/// -[SkipLocalsInit] -public class TrimaVector -{ - private readonly SmaVector _sma1; - private readonly int _count; - private readonly TValue[] _values; - - // Internal state for second stage - private readonly RingBuffer[] _buffers2; - private readonly RingBuffer[] _p_buffers2; - private readonly double[] _lastValidValues2; - - /// - /// Current TRIMA values for all periods. - /// - public ReadOnlySpan Values => _values; - - /// - /// Initializes TrimaVector with specified periods. - /// - /// Array of periods (each must be > 0) - public TrimaVector(int[] periods) - { - _count = periods.Length; - _values = new TValue[_count]; - _buffers2 = new RingBuffer[_count]; - _p_buffers2 = new RingBuffer[_count]; - _lastValidValues2 = new double[_count]; - - int[] p1 = new int[_count]; - - for (int i = 0; i < _count; i++) - { - ArgumentOutOfRangeException.ThrowIfLessThanOrEqual(periods[i], 0); - p1[i] = periods[i] / 2 + 1; - int p2 = (periods[i] + 1) / 2; - - _buffers2[i] = new RingBuffer(p2); - _p_buffers2[i] = new RingBuffer(p2); - } - - _sma1 = new SmaVector(p1); - } - - /// - /// Resets all TRIMA states. - /// - public void Reset() - { - _sma1.Reset(); - for (int i = 0; i < _count; i++) - { - _buffers2[i].Clear(); - _p_buffers2[i].Clear(); - } - Array.Clear(_lastValidValues2); - Array.Clear(_values); - } - - /// - /// Updates TRIMAs with the given value. - /// Uses last-value substitution: invalid inputs (NaN/Infinity) are replaced with - /// the last known good value, providing continuity in the output series. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Array of TRIMA values - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public TValue[] Update(TValue input, bool isNew = true) - { - // First pass: SMA1 - var sma1Results = _sma1.Update(input, isNew); - - // Second pass: SMA2 (TRIMA) - // We need to feed each SMA1 result into the corresponding SMA2 - // Since SmaVector.Update takes a single input, we can't use it directly for vector-to-vector - // However, SmaVector is designed for single input -> multiple periods - // Here we have multiple inputs (from SMA1) -> multiple periods (for SMA2) - // This means we need to update each SMA2 individually, but SmaVector doesn't support that directly - // Wait, SmaVector structure is: one input -> N periods. - // Here we have N inputs (one for each period from SMA1) -> N periods (one for each period in SMA2). - // So we can't use a single SmaVector for the second stage if the inputs are different. - // We need N separate SMAs for the second stage, OR we need to modify SmaVector to support vector input. - // But wait, TrimaVector is supposed to be optimized. - // Let's look at how we can implement this efficiently. - - // Actually, since each period in TRIMA maps to a specific pair of (p1, p2), - // and the input to the second SMA depends on the output of the first SMA, - // the inputs to the second stage are indeed all different. - // So we can't use SmaVector for the second stage in the same way (single input broadcast to all). - - // We have two options: - // 1. Use an array of Sma objects for the second stage. - // 2. Implement a custom vector-input SMA logic here. - - // Given the goal of high performance and vectorization, option 2 is better but more complex. - // However, for now, to match the structure and ensure correctness, let's use the fact that - // we already have SmaVector which is optimized for ring buffers. - // But SmaVector assumes a single input value for all buffers. - // Here, _sma1 produces an array of values, one for each period. - // _sma2 needs to take these DIFFERENT values. - - // So, we cannot use SmaVector for the second stage if it only supports single input. - // Let's check SmaVector again. Yes, Update takes `TValue input`. - - // So we need to implement the second stage manually using RingBuffers, similar to SmaVector - // but accepting a vector of inputs. - - // Let's refactor: - // Instead of using _sma2 as SmaVector, we'll manage the second stage buffers directly here. - // This duplicates some logic from SmaVector but allows vector-to-vector processing. - - // Actually, since we are implementing TrimaVector, maybe we should just use arrays of RingBuffers - // for both stages directly, to avoid the mismatch. - // But _sma1 is fine because it takes the single external input. - // It's only the second stage that is problematic. - - // Let's implement the second stage buffers directly. - - // Wait, I can't change the class structure mid-method. - // I will implement the class using _sma1 for the first stage, and manual buffers for the second stage. - - // Re-reading my own thought process: - // _sma1.Update(input) returns TValue[] with results for each period. - // We need to feed result[i] into buffer2[i]. - - return UpdateInternal(sma1Results, isNew); - } - - private TValue[] UpdateInternal(TValue[] inputs, bool isNew) - { - if (isNew) - { - for (int i = 0; i < _count; i++) - { - _p_buffers2[i].CopyFrom(_buffers2[i]); - } - } - else - { - for (int i = 0; i < _count; i++) - { - _buffers2[i].CopyFrom(_p_buffers2[i]); - } - } - - for (int i = 0; i < _count; i++) - { - double val = inputs[i].Value; - - // Last-value substitution for the second stage - if (double.IsFinite(val)) - { - _lastValidValues2[i] = val; - } - else - { - val = _lastValidValues2[i]; - } - - _buffers2[i].Add(val); - _values[i] = new TValue(inputs[i].Time, _buffers2[i].Average); - } - - return _values; - } - - /// - /// Calculates TRIMAs for the entire series. - /// - /// Input series - /// Array of TRIMA series - public TSeries[] Calculate(TSeries source) - { - // We can use the Update method for simplicity and correctness, - // or implement a batch calculation for performance. - // Given the complexity of double smoothing, using Update in a loop is safer and cleaner. - // SmaVector.Calculate is optimized, but we have the two-stage issue. - - // Let's use the Update loop approach for now to ensure correctness. - // It will be reasonably fast. - - int len = source.Count; - var resultSeries = new TSeries[_count]; - - // Pre-allocate lists - var tLists = new List[_count]; - var vLists = new List[_count]; - - for (int i = 0; i < _count; i++) - { - tLists[i] = new List(len); - vLists[i] = new List(len); - CollectionsMarshal.SetCount(tLists[i], len); - CollectionsMarshal.SetCount(vLists[i], len); - } - - Reset(); - - for (int t = 0; t < len; t++) - { - var tVal = new TValue(source.Times[t], source.Values[t]); - var results = Update(tVal, isNew: true); - - for (int i = 0; i < _count; i++) - { - CollectionsMarshal.AsSpan(tLists[i])[t] = results[i].Time; - CollectionsMarshal.AsSpan(vLists[i])[t] = results[i].Value; - } - } - - for (int i = 0; i < _count; i++) - { - resultSeries[i] = new TSeries(tLists[i], vLists[i]); - } - - return resultSeries; - } - - /// - /// Calculates TRIMAs for the entire series using specified periods. - /// - /// Input series - /// Array of periods - /// Array of TRIMA series - public static TSeries[] Calculate(TSeries source, int[] periods) - { - var trimaVector = new TrimaVector(periods); - return trimaVector.Calculate(source); - } -} diff --git a/lib/averages/wma/Wma.Notebook.dib b/lib/averages/wma/Wma.Notebook.dib index 9af05684..c6ad308f 100644 --- a/lib/averages/wma/Wma.Notebook.dib +++ b/lib/averages/wma/Wma.Notebook.dib @@ -21,7 +21,7 @@ This notebook demonstrates: 1. **Manual Data Processing**: Understanding Batch vs. Streaming modes. 2. **Streaming with `isNew`**: Handling intra-bar updates. 3. **Large Dataset Processing**: Using Geometric Brownian Motion (GBM) generated data. -4. **Vectorized Operations**: Calculating multiple WMAs simultaneously. +4. **Handling Invalid Values**: Last-value substitution for NaN/Infinity. 5. **WMA vs SMA vs EMA**: Comparing different moving averages. #!csharp @@ -190,63 +190,9 @@ Console.WriteLine($"Match: {Math.Abs(batchLargeResult.Last().Value - lastStreamV #!markdown -## 4. Vectorized WMA (Multiple Periods) +## 4. Handling Invalid Values (NaN/Infinity) -`WmaVector` allows calculating multiple WMAs (e.g., 5, 10, 20) simultaneously. This is useful for comparing different timeframes. - -### Vectorized Batch - -#!csharp - -int[] periods = { 5, 10, 20 }; -Console.WriteLine($"\n--- Vectorized Batch WMA (Periods: {string.Join(", ", periods)}) ---"); - -var wmaVectorBatch = new WmaVector(periods); -var vectorBatchResults = wmaVectorBatch.Calculate(closeSeries); - -for (int i = 0; i < periods.Length; i++) -{ - Console.WriteLine($"WMA({periods[i]}) Last Value: {vectorBatchResults[i].Last().Value:F4}"); -} - -#!markdown - -### Vectorized Streaming - -#!csharp - -Console.WriteLine($"\n--- Vectorized Streaming WMA (Periods: {string.Join(", ", periods)}) ---"); - -var wmaVectorStream = new WmaVector(periods); -TValue[] lastVectorVal = null; - -foreach(var item in closeSeries) -{ - lastVectorVal = wmaVectorStream.Update(item); -} - -for (int i = 0; i < periods.Length; i++) -{ - Console.WriteLine($"WMA({periods[i]}) Last Value: {lastVectorVal[i].Value:F4}"); -} - -// Verification -bool allMatch = true; -for (int i = 0; i < periods.Length; i++) -{ - if (Math.Abs(vectorBatchResults[i].Last().Value - lastVectorVal[i].Value) > 1e-10) - { - allMatch = false; - break; - } -} -Console.WriteLine($"\nAll Vectorized Stream/Batch values match: {allMatch}"); - -#!markdown - -## 5. Handling Invalid Values (NaN/Infinity) - -Both `Wma` and `WmaVector` use **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. +`Wma` uses **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. #!csharp @@ -299,7 +245,7 @@ for (int i = 0; i < seriesWithNaN.Count; i++) #!markdown -## 6. WMA vs SMA vs EMA Comparison +## 5. WMA vs SMA vs EMA Comparison The WMA, SMA, and EMA are all trend-following indicators, but they weight data differently: @@ -346,7 +292,7 @@ Console.WriteLine("- WMA provides a balance between SMA's stability and EMA's re #!markdown -## 7. WMA Weights More Recent Values +## 6. WMA Weights More Recent Values This example demonstrates how WMA weights more recent values compared to SMA. diff --git a/lib/averages/wma/Wma.cs b/lib/averages/wma/Wma.cs index 5e83ad5a..f9c69a60 100644 --- a/lib/averages/wma/Wma.cs +++ b/lib/averages/wma/Wma.cs @@ -13,29 +13,12 @@ namespace QuanTAlib; /// WMA applies linear weighting to data points, giving more weight to recent values. /// Uses dual running sums for O(1) complexity per update. /// -/// Key characteristics: -/// - Linear weighting: newest value has weight n, oldest has weight 1 -/// - More responsive than SMA due to emphasis on recent data -/// - Less lag than SMA, but more than EMA -/// - O(1) time complexity for both update and bar correction -/// - O(1) space complexity for state save/restore (scalars only) +/// Calculation: +/// WMA = (n*P_n + (n-1)*P_(n-1) + ... + 1*P_1) / (n*(n+1)/2) /// -/// Calculation method: -/// WMA = (n*P_n + (n-1)*P_(n-1) + ... + 2*P_2 + 1*P_1) / (n*(n+1)/2) -/// -/// O(1) update formula: +/// O(1) update: /// S_new = S - oldest + newest /// W_new = W - S_old + n*newest -/// WMA = W_new / divisor -/// -/// Bar correction (isNew=false): -/// - Restores to state after last isNew=true -/// - Then replaces the last value with new correction value -/// - All O(1) using scalar state -/// -/// Sources: -/// - https://www.investopedia.com/terms/w/weightedaverage.asp -/// - https://school.stockcharts.com/doku.php?id=technical_indicators:weighted_moving_average /// [SkipLocalsInit] public sealed class Wma @@ -44,32 +27,19 @@ public sealed class Wma private readonly double _divisor; private readonly RingBuffer _buffer; - // Dual running sums for O(1) WMA calculation - private double _sum; // Simple sum of values in window - private double _wsum; // Weighted sum of values in window - private double _p_sum; // Sum AFTER last isNew=true (for correction restore) - private double _p_wsum; // Weighted sum AFTER last isNew=true - private double _p_lastInput; // Input that was added on last isNew=true - private double _lastValidValue; - private double _p_lastValidValue; - private int _tickCount; // Counter for periodic sum resync - - // Resync interval: recalculate sum from buffer every N ticks to prevent drift + private double _sum, _wsum; + private double _p_sum, _p_wsum, _p_lastInput; + private double _lastValidValue, _p_lastValidValue; + private int _tickCount; private const int ResyncInterval = 1000; - /// - /// Display name for the indicator. - /// public string Name { get; } + public TValue Value { get; private set; } + public bool IsHot => _buffer.IsFull; - /// - /// Creates WMA with specified period. - /// - /// Number of values to average (must be > 0) public Wma(int period) { - if (period <= 0) - throw new ArgumentException("Period must be greater than 0", nameof(period)); + if (period <= 0) throw new ArgumentException("Period must be greater than 0", nameof(period)); _period = period; _divisor = period * (period + 1) * 0.5; @@ -77,20 +47,6 @@ public sealed class Wma Name = $"Wma({period})"; } - /// - /// Current WMA value. - /// - public TValue Value { get; private set; } - - /// - /// True if the WMA has enough data to produce valid results. - /// WMA is "hot" when the buffer is full (has received at least 'period' values). - /// - public bool IsHot => _buffer.IsFull; - - /// - /// Gets a valid input value, using last-value substitution for non-finite inputs. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private double GetValidValue(double input) { @@ -102,33 +58,25 @@ public sealed class Wma return _lastValidValue; } - /// - /// Updates internal state with a new value. - /// Shared logic for both streaming and batch-reconstruction. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private void UpdateState(double val) { if (_buffer.IsFull) { - // Buffer is full: O(1) update using dual running sums - double oldSum = _sum; // Capture before update + double oldSum = _sum; double oldest = _buffer.Oldest; _sum = _sum - oldest + val; _wsum = _wsum - oldSum + (_period * val); } else { - // Warmup phase: incrementally build sums int count = _buffer.Count + 1; _sum += val; _wsum += count * val; } - // Update buffer _buffer.Add(val); - // Periodic resync: recalculate sums from scratch to eliminate floating-point drift _tickCount++; if (_buffer.IsFull && _tickCount >= ResyncInterval) { @@ -147,24 +95,14 @@ public sealed class Wma } } - /// - /// Updates WMA with the given value. - /// O(1) for both isNew=true and isNew=false. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Current WMA value [MethodImpl(MethodImplOptions.AggressiveInlining)] public TValue Update(TValue input, bool isNew = true) { if (isNew) { - // Get valid value (this may update _lastValidValue) double val = GetValidValue(input.Value); - UpdateState(val); - // Save state AFTER this update for potential future corrections _p_sum = _sum; _p_wsum = _wsum; _p_lastInput = val; @@ -172,40 +110,24 @@ public sealed class Wma } else { - // Bar correction: restore to state AFTER last isNew=true, then swap last value - // Restore _lastValidValue BEFORE calling GetValidValue _lastValidValue = _p_lastValidValue; - - // Get valid value (this may update _lastValidValue) double val = GetValidValue(input.Value); - // Restore sums to state after last isNew=true _sum = _p_sum; _wsum = _p_wsum; - // Correction: replace _p_lastInput with val - // S_corrected = S - lastInput + val - // W_corrected = W + weight*(val - lastInput), where weight = period (if full) or count (if warmup) int weight = _buffer.IsFull ? _period : _buffer.Count; _sum = _sum - _p_lastInput + val; _wsum += weight * (val - _p_lastInput); - // Update buffer's newest value _buffer.UpdateNewest(val); } - // Calculate WMA using current divisor (handles warmup) double currentDivisor = _buffer.IsFull ? _divisor : _buffer.Count * (_buffer.Count + 1) * 0.5; - double result = _wsum / currentDivisor; - Value = new TValue(input.Time, result); + Value = new TValue(input.Time, _wsum / currentDivisor); return Value; } - /// - /// Updates WMA with the entire series. - /// - /// Input series - /// WMA series public TSeries Update(TSeries source) { if (source.Count == 0) return new TSeries(new List(), new List()); @@ -218,41 +140,30 @@ public sealed class Wma var tSpan = CollectionsMarshal.AsSpan(t); var vSpan = CollectionsMarshal.AsSpan(v); - var sourceValues = source.Values; - var sourceTimes = source.Times; - - // 1. Fast Batch Calculation (SIMD optimized) - Calculate(sourceValues, vSpan, _period); - - // 2. Copy Times - sourceTimes.CopyTo(tSpan); - - // 3. Reconstruct State for subsequent updates - // We need to restore _buffer, _sum, _wsum, and _lastValidValue - // Find the last valid value before the reconstruction window + Calculate(source.Values, vSpan, _period); + source.Times.CopyTo(tSpan); + + // Restore state int windowSize = Math.Min(len, _period); int startIndex = len - windowSize; - // Restore _lastValidValue from before the window if (startIndex > 0) { - // Scan backwards to find last valid value for (int i = startIndex - 1; i >= 0; i--) { - if (double.IsFinite(sourceValues[i])) + if (double.IsFinite(source.Values[i])) { - _lastValidValue = sourceValues[i]; + _lastValidValue = source.Values[i]; break; } } } else { - _lastValidValue = 0; // Reset if starting from 0 + _lastValidValue = 0; } - // Rebuild buffer and sums from last 'period' values using shared logic _buffer.Clear(); _sum = 0; _wsum = 0; @@ -260,41 +171,25 @@ public sealed class Wma for (int i = startIndex; i < len; i++) { - double val = GetValidValue(sourceValues[i]); + double val = GetValidValue(source.Values[i]); UpdateState(val); } - // Save state for potential future corrections _p_sum = _sum; _p_wsum = _wsum; - _p_lastInput = sourceValues[len - 1]; + _p_lastInput = source.Values[len - 1]; _p_lastValidValue = _lastValidValue; Value = new TValue(tSpan[len - 1], vSpan[len - 1]); return new TSeries(t, v); } - /// - /// Calculates WMA for the entire series using a new instance. - /// - /// Input series - /// WMA period - /// WMA series public static TSeries Calculate(TSeries source, int period) { var wma = new Wma(period); return wma.Update(source); } - /// - /// Calculates WMA in-place, writing results to pre-allocated output span. - /// Zero-allocation method for maximum performance. - /// Uses O(1) dual running sum algorithm. - /// Automatically uses SIMD acceleration for large, clean datasets. - /// - /// Input values - /// Output span (must be same length as source) - /// WMA period (must be > 0) [MethodImpl(MethodImplOptions.AggressiveInlining)] public static void Calculate(ReadOnlySpan source, Span output, int period) { @@ -306,8 +201,6 @@ public sealed class Wma int len = source.Length; if (len == 0) return; - // Try SIMD path for large, clean datasets - // Requirements: AVX2 support, large enough dataset, no NaN values const int SimdThreshold = 256; if (Avx2.IsSupported && len >= SimdThreshold && !HasNonFiniteValues(source)) { @@ -318,11 +211,6 @@ public sealed class Wma CalculateScalarCore(source, output, period); } - /// - /// Scalar implementation with NaN handling via last-value substitution. - /// Uses circular buffer for sliding window calculation. - /// Optimized with split loops and periodic resync. - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private static void CalculateScalarCore(ReadOnlySpan source, Span output, int period) { @@ -332,12 +220,10 @@ public sealed class Wma double wsum = 0; double lastValid = 0; - // Ring buffer simulation Span buffer = period <= 512 ? stackalloc double[period] : new double[period]; int bufferIdx = 0; int i = 0; - // Phase 1: Warmup (0 to period-1) int warmupEnd = Math.Min(period, len); for (; i < warmupEnd; i++) { @@ -355,7 +241,6 @@ public sealed class Wma output[i] = wsum / currentDivisor; } - // Phase 2: Hot loop (period to len) int tickCount = 0; for (; i < len; i++) { @@ -365,7 +250,6 @@ public sealed class Wma else val = lastValid; - // O(1) update using dual running sums double oldSum = sum; double oldest = buffer[bufferIdx]; sum = sum - oldest + val; @@ -378,28 +262,17 @@ public sealed class Wma output[i] = wsum / divisor; - // Periodic resync every 1000 ticks tickCount++; if (tickCount >= ResyncInterval) { tickCount = 0; - // Recalculate sums from buffer to prevent drift double recalcSum = 0; double recalcWsum = 0; - // Buffer contains values in order: [oldest ... newest] relative to current bufferIdx - // Actually buffer is circular. - // Oldest is at bufferIdx (which we just wrote to, so it's actually newest now? No, we incremented bufferIdx) - // bufferIdx points to the *next* overwrite location, which holds the *oldest* value. - // So buffer[bufferIdx] is oldest (weight 1). - // buffer[bufferIdx+1] is 2nd oldest (weight 2). - // ... - // buffer[bufferIdx-1] is newest (weight period). for (int k = 0; k < period; k++) { - int idx = (bufferIdx + k) % period; // Use modulo here for simplicity in resync (rare) - // Wait, modulo is slow. - if (idx >= period) idx -= period; // Manual modulo + int idx = bufferIdx + k; + if (idx >= period) idx -= period; double v = buffer[idx]; recalcSum += v; @@ -411,10 +284,6 @@ public sealed class Wma } } - /// - /// SIMD-optimized implementation for WMA calculation. - /// Uses double prefix-sum approach to vectorize the coupled recurrence. - /// [MethodImpl(MethodImplOptions.AggressiveOptimization)] private static unsafe void CalculateSimdCore(ReadOnlySpan source, Span output, int period) { @@ -427,7 +296,6 @@ public sealed class Wma double divisor = period * (period + 1) * 0.5; double invDivisor = 1.0 / divisor; - // Phase 1: Warmup - scalar int warmupEnd = Math.Min(period, len); double sum = 0; double wsum = 0; @@ -443,13 +311,11 @@ public sealed class Wma if (len <= period) return; - // Phase 2: SIMD hot loop var vInvDivisor = Vector256.Create(invDivisor); var vPeriod = Vector256.Create((double)period); var vZero = Vector256.Zero; int simdEnd = period + ((len - period) / VectorWidth) * VectorWidth; - // Initialize vector state var vSumState = Vector256.Create(sum); var vWsumState = Vector256.Create(wsum); @@ -458,26 +324,17 @@ public sealed class Wma { int nextSync = Math.Min(simdEnd, idx + ResyncInterval); - // Inner hot loop without branches - // Unrolled 2x (process 8 doubles per iteration) - // Optimized Parallel Execution: - // - Parallel prefix sums for DeltaS - // - Fast S_shifted calculation using (S - DeltaS) - // - Parallel prefix sums for U int unrolledSync = nextSync - (2 * VectorWidth); for (; idx <= unrolledSync; idx += 2 * VectorWidth) { - // Load data for both iterations var vNew1 = Avx.LoadVector256(srcPtr + idx); var vOld1 = Avx.LoadVector256(srcPtr + idx - period); var vNew2 = Avx.LoadVector256(srcPtr + idx + VectorWidth); var vOld2 = Avx.LoadVector256(srcPtr + idx + VectorWidth - period); - // 1. Update Sum (S) - Parallel Prefix Sums var vDeltaS1 = Avx.Subtract(vNew1, vOld1); var vDeltaS2 = Avx.Subtract(vNew2, vOld2); - // Prefix Sum DeltaS1 var vShiftS1_1 = Avx2.Permute4x64(vDeltaS1.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShiftS1_1 = Avx.Blend(vZero, vShiftS1_1, 0b_1110); var vPS_DeltaS1 = Avx.Add(vDeltaS1, vShiftS1_1); @@ -485,7 +342,6 @@ public sealed class Wma vShiftS2_1 = Avx.Blend(vZero, vShiftS2_1, 0b_1100); vPS_DeltaS1 = Avx.Add(vPS_DeltaS1, vShiftS2_1); - // Prefix Sum DeltaS2 var vShiftS1_2 = Avx2.Permute4x64(vDeltaS2.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShiftS1_2 = Avx.Blend(vZero, vShiftS1_2, 0b_1110); var vPS_DeltaS2 = Avx.Add(vDeltaS2, vShiftS1_2); @@ -493,14 +349,10 @@ public sealed class Wma vShiftS2_2 = Avx.Blend(vZero, vShiftS2_2, 0b_1100); vPS_DeltaS2 = Avx.Add(vPS_DeltaS2, vShiftS2_2); - // Combine Sums var vSums1 = Avx.Add(vSumState, vPS_DeltaS1); var vLastS1 = Avx2.Permute4x64(vSums1.AsUInt64(), 0b_11_11_11_11).AsDouble(); var vSums2 = Avx.Add(vLastS1, vPS_DeltaS2); - // 2. Update Weighted Sum (W) - // Optimization: S_shifted = S - DeltaS - // This avoids expensive Permute/Blend operations var vSumsShifted1 = Avx.Subtract(vSums1, vDeltaS1); var vSumsShifted2 = Avx.Subtract(vSums2, vDeltaS2); @@ -516,7 +368,6 @@ public sealed class Wma vU2 = Avx.Subtract(Avx.Multiply(vPeriod, vNew2), vSumsShifted2); } - // Prefix Sum W1 var vShiftW1_1 = Avx2.Permute4x64(vU1.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShiftW1_1 = Avx.Blend(vZero, vShiftW1_1, 0b_1110); var vPW1_1 = Avx.Add(vU1, vShiftW1_1); @@ -524,7 +375,6 @@ public sealed class Wma vShiftW2_1 = Avx.Blend(vZero, vShiftW2_1, 0b_1100); var vPW2_1 = Avx.Add(vPW1_1, vShiftW2_1); - // Prefix Sum W2 var vShiftW1_2 = Avx2.Permute4x64(vU2.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShiftW1_2 = Avx.Blend(vZero, vShiftW1_2, 0b_1110); var vPW1_2 = Avx.Add(vU2, vShiftW1_2); @@ -532,32 +382,24 @@ public sealed class Wma vShiftW2_2 = Avx.Blend(vZero, vShiftW2_2, 0b_1100); var vPW2_2 = Avx.Add(vPW1_2, vShiftW2_2); - // Combine Weighted Sums var vWsums1 = Avx.Add(vWsumState, vPW2_1); var vLastW1 = Avx2.Permute4x64(vWsums1.AsUInt64(), 0b_11_11_11_11).AsDouble(); var vWsums2 = Avx.Add(vLastW1, vPW2_2); - // Store results Avx.Store(outPtr + idx, Avx.Multiply(vWsums1, vInvDivisor)); Avx.Store(outPtr + idx + VectorWidth, Avx.Multiply(vWsums2, vInvDivisor)); - // Update state for next iteration vSumState = Avx2.Permute4x64(vSums2.AsUInt64(), 0b_11_11_11_11).AsDouble(); vWsumState = Avx2.Permute4x64(vWsums2.AsUInt64(), 0b_11_11_11_11).AsDouble(); } - // Handle remaining vectors (if any) for (; idx < nextSync; idx += VectorWidth) { - // Load 4 entering values and 4 leaving values var vNew = Avx.LoadVector256(srcPtr + idx); var vOld = Avx.LoadVector256(srcPtr + idx - period); - // 1. Update Sum (S) - // Delta S = New - Old var vDeltaS = Avx.Subtract(vNew, vOld); - // Prefix sum of Delta S var vShiftS1 = Avx2.Permute4x64(vDeltaS.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShiftS1 = Avx.Blend(vZero, vShiftS1, 0b_1110); var vPS1 = Avx.Add(vDeltaS, vShiftS1); @@ -566,19 +408,14 @@ public sealed class Wma vShiftS2 = Avx.Blend(vZero, vShiftS2, 0b_1100); var vPS2 = Avx.Add(vPS1, vShiftS2); - // Add previous sum state var vSums = Avx.Add(vSumState, vPS2); - // 2. Update Weighted Sum (W) - // Shift vSums right and insert sum (S_t) at pos 0 var vSumsShifted = Avx2.Permute4x64(vSums.AsUInt64(), 0b_10_01_00_00).AsDouble(); vSumsShifted = Avx.Blend(vSumState, vSumsShifted, 0b_1110); - // U = (n * New) - S_shifted var vTerm1 = Avx.Multiply(vPeriod, vNew); var vU = Avx.Subtract(vTerm1, vSumsShifted); - // Prefix sum of U var vShiftW1 = Avx2.Permute4x64(vU.AsUInt64(), 0b_10_01_00_00).AsDouble(); vShiftW1 = Avx.Blend(vZero, vShiftW1, 0b_1110); var vPW1 = Avx.Add(vU, vShiftW1); @@ -587,26 +424,17 @@ public sealed class Wma vShiftW2 = Avx.Blend(vZero, vShiftW2, 0b_1100); var vPW2 = Avx.Add(vPW1, vShiftW2); - // Add previous wsum state var vWsums = Avx.Add(vWsumState, vPW2); - // Store result var vResult = Avx.Multiply(vWsums, vInvDivisor); Avx.Store(outPtr + idx, vResult); - // Update state for next iteration vSumState = Avx2.Permute4x64(vSums.AsUInt64(), 0b_11_11_11_11).AsDouble(); vWsumState = Avx2.Permute4x64(vWsums.AsUInt64(), 0b_11_11_11_11).AsDouble(); } - // Periodic resync if (idx < len) { - // Extract scalar state for resync logic - sum = vSumState.GetElement(0); - wsum = vWsumState.GetElement(0); - - // Recalculate sums from scratch int lastIdx = idx - 1; double recalcSum = 0; double recalcWsum = 0; @@ -619,17 +447,14 @@ public sealed class Wma sum = recalcSum; wsum = recalcWsum; - // Update vector state after resync vSumState = Vector256.Create(sum); vWsumState = Vector256.Create(wsum); } } - // Extract final scalar state for tail sum = vSumState.GetElement(0); wsum = vWsumState.GetElement(0); - // Phase 3: Scalar tail for (; idx < len; idx++) { double val = srcPtr[idx]; @@ -642,9 +467,6 @@ public sealed class Wma } } - /// - /// Checks if span contains any non-finite values (NaN or Infinity). - /// [MethodImpl(MethodImplOptions.AggressiveInlining)] private static bool HasNonFiniteValues(ReadOnlySpan span) { @@ -656,19 +478,10 @@ public sealed class Wma return false; } - /// - /// Resets the WMA state. - /// public void Reset() { _buffer.Clear(); - _sum = 0; - _wsum = 0; - _p_sum = 0; - _p_wsum = 0; - _p_lastInput = 0; - _lastValidValue = 0; - _p_lastValidValue = 0; + _sum = _wsum = _p_sum = _p_wsum = _p_lastInput = _lastValidValue = _p_lastValidValue = 0; Value = default; } } diff --git a/lib/averages/wma/Wma.md b/lib/averages/wma/Wma.md index 11aa9483..1f2d3a3c 100644 --- a/lib/averages/wma/Wma.md +++ b/lib/averages/wma/Wma.md @@ -114,33 +114,9 @@ Console.WriteLine($"Last WMA: {wmaOutput[^1]}"); * **O(1) per-bar** via dual running sums * **Compatible** with `ArrayPool` for buffer management -### Multi-Period WMA (`WmaVector`) - -The `WmaVector` class calculates multiple WMAs with different periods on the same input series simultaneously. - -```csharp -using QuanTAlib; - -// Initialize with multiple periods -int[] periods = { 5, 10, 20 }; -var wmaVector = new WmaVector(periods); - -// Streaming update -TValue[] results = wmaVector.Update(new TValue(time, price)); - -// Access values -Console.WriteLine($"WMA(5): {results[0].Value}"); -Console.WriteLine($"WMA(10): {results[1].Value}"); -Console.WriteLine($"WMA(20): {results[2].Value}"); - -// Batch calculation -TSeries source = ...; -TSeries[] seriesResults = wmaVector.Calculate(source); -``` - ### Bar Correction (isNew Parameter) -Both `Wma` and `WmaVector` support intra-bar updates for real-time trading systems: +`Wma` supports intra-bar updates for real-time trading systems: ```csharp var wma = new Wma(10); @@ -166,7 +142,7 @@ wma.Update(new TValue(time + 1, 101.2), isNew: true); ### Handling Invalid Values (NaN/Infinity) -Both `Wma` and `WmaVector` use **last-value substitution** for handling invalid inputs: +`Wma` uses **last-value substitution** for handling invalid inputs: ```csharp var wma = new Wma(10); diff --git a/lib/averages/wma/WmaVector.Tests.cs b/lib/averages/wma/WmaVector.Tests.cs deleted file mode 100644 index 2805001d..00000000 --- a/lib/averages/wma/WmaVector.Tests.cs +++ /dev/null @@ -1,407 +0,0 @@ - -namespace QuanTAlib.Tests; - -public class WmaVectorTests -{ - [Fact] - public void Initialization_WithPeriods_Works() - { - int[] periods = { 5, 10, 20 }; - var wmaVector = new WmaVector(periods); - - var res = wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(3, res.Length); - Assert.Equal(100.0, res[0].Value, 1e-9); - Assert.Equal(100.0, res[1].Value, 1e-9); - Assert.Equal(100.0, res[2].Value, 1e-9); - } - - [Fact] - public void Initialization_WithZeroPeriod_ThrowsArgumentException() - { - int[] periods = { 10, 0, 20 }; - - Assert.Throws(() => new WmaVector(periods)); - } - - [Fact] - public void Initialization_WithNegativePeriod_ThrowsArgumentException() - { - int[] periods = { 10, -5, 20 }; - - Assert.Throws(() => new WmaVector(periods)); - } - - [Fact] - public void Calc_Streaming_MatchesSingleWma() - { - int[] periods = { 5, 10, 20 }; - var wmaVector = new WmaVector(periods); - var wmaSingles = periods.Select(p => new Wma(p)).ToArray(); - - var values = new double[] { 10, 20, 30, 40, 50, 40, 30, 20, 10 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = wmaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = wmaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - Assert.Equal(singleRes.Time, multiRes[i].Time); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Calc_Series_MatchesSingleWma() - { - int[] periods = { 5, 10, 20 }; - var wmaVector = new WmaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var multiRes = wmaVector.Calculate(series); - - // Reset and recalculate for comparison - var wmaSingles = periods.Select(p => new Wma(p)).ToArray(); - for (int j = 0; j < len; j++) - { - var tVal = new TValue(new DateTime(t[j], DateTimeKind.Utc), v[j]); - for (int i = 0; i < periods.Length; i++) - { - var singleRes = wmaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Values[j], 1e-8); - } - } - } - - [Fact] - public void Calc_Series_MatchesStreaming() - { - int[] periods = { 5, 10, 20 }; - var wmaVectorBatch = new WmaVector(periods); - var wmaVectorStream = new WmaVector(periods); - - int len = 100; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var batchRes = wmaVectorBatch.Calculate(series); - - for (int i = 0; i < len; i++) - { - var tVal = new TValue(new DateTime(t[i], DateTimeKind.Utc), v[i]); - var streamRes = wmaVectorStream.Update(tVal); - - for (int j = 0; j < periods.Length; j++) - { - Assert.Equal(batchRes[j].Values[i], streamRes[j].Value, 1e-9); - } - } - } - - [Fact] - public void Calculate_Static_MatchesInstanceMethod() - { - int[] periods = { 5, 10, 20 }; - - int len = 50; - var t = new System.Collections.Generic.List(len); - var v = new System.Collections.Generic.List(len); - var now = DateTime.UtcNow; - - for (int i = 0; i < len; i++) - { - t.Add(now.AddMinutes(i).Ticks); - v.Add(Math.Sin(i * 0.1) * 100); - } - - var series = new TSeries(t, v); - - var instanceWma = new WmaVector(periods); - var instanceRes = instanceWma.Calculate(series); - - var staticRes = WmaVector.Calculate(series, periods); - - for (int i = 0; i < periods.Length; i++) - { - Assert.Equal(instanceRes[i].Count, staticRes[i].Count); - for (int j = 0; j < len; j++) - { - Assert.Equal(instanceRes[i].Values[j], staticRes[i].Values[j], 1e-9); - } - } - } - - [Fact] - public void Reset_ClearsState() - { - int[] periods = { 10 }; - var wmaVector = new WmaVector(periods); - - wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - wmaVector.Update(new TValue(DateTime.UtcNow, 200.0)); - wmaVector.Reset(); - - var res = wmaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - - Assert.Equal(50.0, res[0].Value, 1e-9); - } - - [Fact] - public void Update_NaN_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var wmaVector = new WmaVector(periods); - - wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - wmaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterNaN = wmaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in resultAfterNaN) - { - Assert.True(double.IsFinite(result.Value), $"Expected finite value but got {result.Value}"); - } - } - - [Fact] - public void Update_Infinity_Input_UsesLastValidValue() - { - int[] periods = { 10, 20 }; - var wmaVector = new WmaVector(periods); - - wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - wmaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - - var resultAfterPosInf = wmaVector.Update(new TValue(DateTime.UtcNow, double.PositiveInfinity)); - foreach (var result in resultAfterPosInf) - { - Assert.True(double.IsFinite(result.Value)); - } - - var resultAfterNegInf = wmaVector.Update(new TValue(DateTime.UtcNow, double.NegativeInfinity)); - foreach (var result in resultAfterNegInf) - { - Assert.True(double.IsFinite(result.Value)); - } - } - - [Fact] - public void Update_MultipleNaN_ContinuesWithLastValid() - { - int[] periods = { 5, 10 }; - var wmaVector = new WmaVector(periods); - - wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - wmaVector.Update(new TValue(DateTime.UtcNow, 110.0)); - wmaVector.Update(new TValue(DateTime.UtcNow, 120.0)); - - var r1 = wmaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r2 = wmaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - var r3 = wmaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - foreach (var result in r1) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r2) Assert.True(double.IsFinite(result.Value)); - foreach (var result in r3) Assert.True(double.IsFinite(result.Value)); - } - - [Fact] - public void Calculate_Series_HandlesNaN() - { - int[] periods = { 5, 10 }; - var wmaVector = new WmaVector(periods); - - var t = new System.Collections.Generic.List(); - var v = new System.Collections.Generic.List(); - var now = DateTime.UtcNow; - - t.Add(now.Ticks); v.Add(100.0); - t.Add(now.AddMinutes(1).Ticks); v.Add(110.0); - t.Add(now.AddMinutes(2).Ticks); v.Add(double.NaN); - t.Add(now.AddMinutes(3).Ticks); v.Add(120.0); - t.Add(now.AddMinutes(4).Ticks); v.Add(double.PositiveInfinity); - t.Add(now.AddMinutes(5).Ticks); v.Add(130.0); - - var series = new TSeries(t, v); - var results = wmaVector.Calculate(series); - - foreach (var periodResults in results) - { - foreach (var val in periodResults.Values) - { - Assert.True(double.IsFinite(val), $"Expected finite value but got {val}"); - } - } - } - - [Fact] - public void Reset_ClearsLastValidValue() - { - int[] periods = { 10 }; - var wmaVector = new WmaVector(periods); - - wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - wmaVector.Update(new TValue(DateTime.UtcNow, double.NaN)); - - wmaVector.Reset(); - - var result = wmaVector.Update(new TValue(DateTime.UtcNow, 50.0)); - Assert.Equal(50.0, result[0].Value, 1e-9); - } - - [Fact] - public void NaN_Handling_MatchesSingleWma() - { - int[] periods = { 5, 10, 20 }; - var wmaVector = new WmaVector(periods); - var wmaSingles = periods.Select(p => new Wma(p)).ToArray(); - - var values = new double[] { 10, 20, double.NaN, 40, double.PositiveInfinity, 60, 70 }; - var time = DateTime.UtcNow; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - var multiRes = wmaVector.Update(tVal); - - for (int i = 0; i < periods.Length; i++) - { - var singleRes = wmaSingles[i].Update(tVal); - Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); - } - - time = time.AddMinutes(1); - } - } - - [Fact] - public void Values_Property_UpdatesAfterUpdate() - { - int[] periods = { 5, 10 }; - var wmaVector = new WmaVector(periods); - - var result = wmaVector.Update(new TValue(DateTime.UtcNow, 100.0)); - - Assert.Equal(result[0].Value, wmaVector.Values[0].Value); - Assert.Equal(result[1].Value, wmaVector.Values[1].Value); - } - - [Fact] - public void Values_Property_UpdatesAfterCalculate() - { - int[] periods = { 5, 10 }; - var wmaVector = new WmaVector(periods); - - var t = new System.Collections.Generic.List { 100, 200, 300 }; - var v = new System.Collections.Generic.List { 10.0, 20.0, 30.0 }; - var series = new TSeries(t, v); - - var results = wmaVector.Calculate(series); - - Assert.Equal(results[0].Last.Value, wmaVector.Values[0].Value, 1e-9); - Assert.Equal(results[1].Last.Value, wmaVector.Values[1].Value, 1e-9); - } - - [Fact] - public void Update_BarCorrection_WorksCorrectly() - { - int[] periods = { 3 }; - var wmaVector = new WmaVector(periods); - - wmaVector.Update(new TValue(DateTime.UtcNow, 10.0), isNew: true); - wmaVector.Update(new TValue(DateTime.UtcNow, 20.0), isNew: true); - wmaVector.Update(new TValue(DateTime.UtcNow, 30.0), isNew: true); - - // WMA(3) of 10,20,30 = (1*10 + 2*20 + 3*30) / 6 = 140/6 = 23.333... - var res1 = wmaVector.Values[0].Value; - Assert.Equal(140.0 / 6.0, res1, 1e-9); - - // Correct the last bar to 60 - var res2 = wmaVector.Update(new TValue(DateTime.UtcNow, 60.0), isNew: false); - - // WMA(3) of 10,20,60 = (1*10 + 2*20 + 3*60) / 6 = (10 + 40 + 180) / 6 = 230/6 = 38.333... - Assert.Equal(230.0 / 6.0, res2[0].Value, 1e-9); - } - - [Fact] - public void WMA_MatchesExpectedValues() - { - int[] periods = { 3 }; - var wmaVector = new WmaVector(periods); - - // Test sequence: 10, 20, 30, 40, 50 - // WMA(3) weights: [1, 2, 3], divisor = 6 - // Bar 1: 10 (only value) = 10 - // Bar 2: (1*10 + 2*20) / 3 = 50/3 = 16.666... - // Bar 3: (1*10 + 2*20 + 3*30) / 6 = 140/6 = 23.333... - // Bar 4: (1*20 + 2*30 + 3*40) / 6 = 200/6 = 33.333... - // Bar 5: (1*30 + 2*40 + 3*50) / 6 = 260/6 = 43.333... - double[] expected = [10.0, 50.0/3.0, 140.0/6.0, 200.0/6.0, 260.0/6.0]; - var values = new double[] { 10, 20, 30, 40, 50 }; - var time = DateTime.UtcNow; - - for (int i = 0; i < values.Length; i++) - { - var res = wmaVector.Update(new TValue(time, values[i])); - Assert.Equal(expected[i], res[0].Value, 1e-9); - time = time.AddMinutes(1); - } - } - - [Fact] - public void WMA_MoreWeightOnRecentValues() - { - int[] periods = { 3 }; - var wmaVector = new WmaVector(periods); - var smaVector = new SmaVector(periods); - - var values = new double[] { 10, 20, 100 }; // High recent value - var time = DateTime.UtcNow; - - TValue[] wmaRes = null!; - TValue[] smaRes = null!; - - foreach (var val in values) - { - var tVal = new TValue(time, val); - wmaRes = wmaVector.Update(tVal); - smaRes = smaVector.Update(tVal); - time = time.AddMinutes(1); - } - - // WMA should be higher than SMA because it weights the high recent value more - // SMA = (10 + 20 + 100) / 3 = 43.333... - // WMA = (1*10 + 2*20 + 3*100) / 6 = (10 + 40 + 300) / 6 = 58.333... - Assert.True(wmaRes[0].Value > smaRes[0].Value); - Assert.Equal(350.0 / 6.0, wmaRes[0].Value, 1e-9); - Assert.Equal(130.0 / 3.0, smaRes[0].Value, 1e-9); - } -} diff --git a/lib/averages/wma/WmaVector.cs b/lib/averages/wma/WmaVector.cs deleted file mode 100644 index 72d8b480..00000000 --- a/lib/averages/wma/WmaVector.cs +++ /dev/null @@ -1,264 +0,0 @@ -using System.Runtime.CompilerServices; -using System.Runtime.InteropServices; - -namespace QuanTAlib; - -/// -/// Multi-Period Weighted Moving Average (WMA) - O(1) optimized per period. -/// Calculates multiple WMAs with different periods for the same input series. -/// Uses dual running sums for O(1) complexity per update per period. -/// Uses last-value substitution for invalid inputs (NaN/Infinity). -/// -[SkipLocalsInit] -public class WmaVector -{ - private readonly int[] _periods; - private readonly double[] _divisors; - private readonly RingBuffer[] _buffers; - private readonly double[] _sums; // Simple sums for each period - private readonly double[] _wsums; // Weighted sums for each period - private readonly double[] _p_sums; // Saved simple sums for bar correction - private readonly double[] _p_wsums; // Saved weighted sums for bar correction - private readonly double[] _p_lastInputs; // Last inputs for bar correction - private readonly int _count; - private double _lastValidValue; - private double _p_lastValidValue; - - /// - /// Current WMA values for all periods. - /// - public ReadOnlySpan Values => _values; - - private readonly TValue[] _values; - - /// - /// Initializes WmaVector with specified periods. - /// - /// Array of periods (each must be > 0) - public WmaVector(int[] periods) - { - _count = periods.Length; - _periods = new int[_count]; - _divisors = new double[_count]; - _buffers = new RingBuffer[_count]; - _sums = new double[_count]; - _wsums = new double[_count]; - _p_sums = new double[_count]; - _p_wsums = new double[_count]; - _p_lastInputs = new double[_count]; - _values = new TValue[_count]; - - for (int i = 0; i < _count; i++) - { - ArgumentOutOfRangeException.ThrowIfLessThanOrEqual(periods[i], 0); - _periods[i] = periods[i]; - _divisors[i] = periods[i] * (periods[i] + 1) * 0.5; - _buffers[i] = new RingBuffer(periods[i]); - } - } - - /// - /// Gets a valid input value, using last-value substitution for non-finite inputs. - /// - [MethodImpl(MethodImplOptions.AggressiveInlining)] - private double GetValidValue(double input) - { - if (double.IsFinite(input)) - { - _lastValidValue = input; - return input; - } - return _lastValidValue; - } - - /// - /// Resets all WMA states. - /// - public void Reset() - { - for (int i = 0; i < _count; i++) - { - _buffers[i].Clear(); - _sums[i] = 0; - _wsums[i] = 0; - _p_sums[i] = 0; - _p_wsums[i] = 0; - _p_lastInputs[i] = 0; - } - _lastValidValue = 0; - _p_lastValidValue = 0; - Array.Clear(_values); - } - - /// - /// Updates WMAs with the given value. - /// Uses last-value substitution: invalid inputs (NaN/Infinity) are replaced with - /// the last known good value, providing continuity in the output series. - /// O(1) complexity per period using dual running sums. - /// - /// Input value - /// True for new bar, false for update to current bar (default: true) - /// Array of WMA values - [MethodImpl(MethodImplOptions.AggressiveInlining)] - public TValue[] Update(TValue input, bool isNew = true) - { - if (isNew) - { - // Get valid value (this may update _lastValidValue) - double val = GetValidValue(input.Value); - - for (int i = 0; i < _count; i++) - { - int period = _periods[i]; - var buffer = _buffers[i]; - - if (buffer.IsFull) - { - // Buffer is full: O(1) update using dual running sums - double oldSum = _sums[i]; - double oldest = buffer.Oldest; - _sums[i] = _sums[i] - oldest + val; - _wsums[i] = _wsums[i] - oldSum + (period * val); - } - else - { - // Warmup phase: incrementally build sums - int count = buffer.Count + 1; - _sums[i] += val; - _wsums[i] += count * val; - } - - buffer.Add(val); - - // Save state AFTER this update for potential future corrections - _p_sums[i] = _sums[i]; - _p_wsums[i] = _wsums[i]; - _p_lastInputs[i] = val; - - // Calculate WMA - double currentDivisor = buffer.IsFull ? _divisors[i] : buffer.Count * (buffer.Count + 1) * 0.5; - _values[i] = new TValue(input.Time, _wsums[i] / currentDivisor); - } - - _p_lastValidValue = _lastValidValue; - } - else - { - // Bar correction: restore to state AFTER last isNew=true, then swap last value - _lastValidValue = _p_lastValidValue; - double val = GetValidValue(input.Value); - - for (int i = 0; i < _count; i++) - { - int period = _periods[i]; - var buffer = _buffers[i]; - - // Restore sums to state after last isNew=true - _sums[i] = _p_sums[i]; - _wsums[i] = _p_wsums[i]; - - // Correction: replace _p_lastInputs[i] with val - int weight = buffer.IsFull ? period : buffer.Count; - _sums[i] = _sums[i] - _p_lastInputs[i] + val; - _wsums[i] += weight * (val - _p_lastInputs[i]); - - buffer.UpdateNewest(val); - - // Calculate WMA - double currentDivisor = buffer.IsFull ? _divisors[i] : buffer.Count * (buffer.Count + 1) * 0.5; - _values[i] = new TValue(input.Time, _wsums[i] / currentDivisor); - } - } - - return _values; - } - - /// - /// Calculates WMAs for the entire series. - /// - /// Input series - /// Array of WMA series - public TSeries[] Calculate(TSeries source) - { - int len = source.Count; - var resultSeries = new TSeries[_count]; - - // Reset state for fresh calculation - Reset(); - - // Pre-allocate lists - var tLists = new List[_count]; - var vLists = new List[_count]; - - for (int i = 0; i < _count; i++) - { - tLists[i] = new List(len); - vLists[i] = new List(len); - CollectionsMarshal.SetCount(tLists[i], len); - CollectionsMarshal.SetCount(vLists[i], len); - } - - var sourceValues = source.Values; - var sourceTimes = source.Times; - - for (int t = 0; t < len; t++) - { - double val = sourceValues[t]; - long time = sourceTimes[t]; - - // Last-value substitution: replace non-finite inputs with last valid value - val = GetValidValue(val); - - for (int i = 0; i < _count; i++) - { - int period = _periods[i]; - var buffer = _buffers[i]; - - if (buffer.IsFull) - { - // Buffer is full: O(1) update - double oldSum = _sums[i]; - double oldest = buffer.Oldest; - _sums[i] = _sums[i] - oldest + val; - _wsums[i] = _wsums[i] - oldSum + (period * val); - } - else - { - // Warmup phase - int count = buffer.Count + 1; - _sums[i] += val; - _wsums[i] += count * val; - } - - buffer.Add(val); - - CollectionsMarshal.AsSpan(tLists[i])[t] = time; - double currentDivisor = buffer.IsFull ? _divisors[i] : buffer.Count * (buffer.Count + 1) * 0.5; - CollectionsMarshal.AsSpan(vLists[i])[t] = _wsums[i] / currentDivisor; - } - } - - // Create TSeries and update Values - for (int i = 0; i < _count; i++) - { - resultSeries[i] = new TSeries(tLists[i], vLists[i]); - var lastT = CollectionsMarshal.AsSpan(tLists[i])[len - 1]; - var lastV = CollectionsMarshal.AsSpan(vLists[i])[len - 1]; - _values[i] = new TValue(lastT, lastV); - } - - return resultSeries; - } - - /// - /// Calculates WMAs for the entire series using specified periods. - /// - /// Input series - /// Array of periods - /// Array of WMA series - public static TSeries[] Calculate(TSeries source, int[] periods) - { - var wmaVector = new WmaVector(periods); - return wmaVector.Calculate(source); - } -} diff --git a/lib/core/ringbuffer/RingBuffer.cs b/lib/core/ringbuffer/RingBuffer.cs index 6a6ba491..9fed9aa3 100644 --- a/lib/core/ringbuffer/RingBuffer.cs +++ b/lib/core/ringbuffer/RingBuffer.cs @@ -23,9 +23,9 @@ public sealed class RingBuffer : IEnumerable { private readonly double[] _buffer; private readonly int _capacity; - private int _head; // Next write position (also start position when full) - private int _count; // Current number of elements - private double _sum; // Running sum of all elements + private int _head; + private int _count; + private double _sum; /// /// Creates a new RingBuffer with the specified capacity. @@ -114,7 +114,6 @@ public sealed class RingBuffer : IEnumerable get { if (_count == 0) return 0; - // When full, _head points to oldest; otherwise start is 0 int start = _count == _capacity ? _head : 0; return _buffer[start]; } @@ -143,7 +142,6 @@ public sealed class RingBuffer : IEnumerable if (_count == _capacity) { - // Buffer is full: remove oldest value from sum removed = _buffer[_head]; _sum -= removed; } @@ -243,13 +241,11 @@ public sealed class RingBuffer : IEnumerable int start = _count == _capacity ? _head : 0; - // Check if contiguous (no wrap) if (start + _count <= _capacity) { return new ReadOnlySpan(_buffer, start, _count); } - // Wrapped - need to copy return new ReadOnlySpan(ToArray()); } diff --git a/lib/core/simd/SimdExtensions.cs b/lib/core/simd/SimdExtensions.cs index 250cc4aa..9cf53d30 100644 --- a/lib/core/simd/SimdExtensions.cs +++ b/lib/core/simd/SimdExtensions.cs @@ -112,7 +112,6 @@ public static class SimdExtensions return true; } - // Check remaining elements with scalar for (; i < span.Length; i++) { if (!double.IsFinite(span[i])) @@ -144,19 +143,16 @@ public static class SimdExtensions int vectorSize = Vector.Count; int i = 0; - // Process in vector chunks for (; i <= span.Length - vectorSize; i += vectorSize) { var vector = new Vector(span.Slice(i, vectorSize)); sum += vector; } - // Horizontal sum of vector double result = 0.0; for (int j = 0; j < vectorSize; j++) result += sum[j]; - // Process remaining elements for (; i < span.Length; i++) result += span[i]; @@ -186,14 +182,12 @@ public static class SimdExtensions var minVec = new Vector(span.Slice(0, vectorSize)); int i = vectorSize; - // Process in vector chunks for (; i <= span.Length - vectorSize; i += vectorSize) { var vector = new Vector(span.Slice(i, vectorSize)); minVec = Vector.Min(minVec, vector); } - // Find minimum within vector double result = minVec[0]; for (int j = 1; j < vectorSize; j++) { @@ -201,7 +195,6 @@ public static class SimdExtensions result = minVec[j]; } - // Process remaining elements for (; i < span.Length; i++) { if (span[i] < result) @@ -234,14 +227,12 @@ public static class SimdExtensions var maxVec = new Vector(span.Slice(0, vectorSize)); int i = vectorSize; - // Process in vector chunks for (; i <= span.Length - vectorSize; i += vectorSize) { var vector = new Vector(span.Slice(i, vectorSize)); maxVec = Vector.Max(maxVec, vector); } - // Find maximum within vector double result = maxVec[0]; for (int j = 1; j < vectorSize; j++) { @@ -249,7 +240,6 @@ public static class SimdExtensions result = maxVec[j]; } - // Process remaining elements for (; i < span.Length; i++) { if (span[i] > result) @@ -285,12 +275,17 @@ public static class SimdExtensions { if (span.Length < 2) return double.NaN; - // Guard against non-finite inputs - if (span.ContainsNonFinite()) return double.NaN; + double m; + if (mean.HasValue) + { + if (span.ContainsNonFinite()) return double.NaN; + m = mean.Value; + } + else + { + m = span.AverageSIMD(); + } - double m = mean ?? span.AverageSIMD(); - - // Guard against non-finite mean (could be passed in or computed from non-finite values) if (!double.IsFinite(m)) return double.NaN; if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) @@ -300,7 +295,6 @@ public static class SimdExtensions int vectorSize = Vector.Count; int i = 0; - // Process in vector chunks for (; i <= span.Length - vectorSize; i += vectorSize) { var vector = new Vector(span.Slice(i, vectorSize)); @@ -308,12 +302,10 @@ public static class SimdExtensions sumSq += diff * diff; } - // Horizontal sum of vector double result = 0.0; for (int j = 0; j < vectorSize; j++) result += sumSq[j]; - // Process remaining elements for (; i < span.Length; i++) { double diff = span[i] - m; @@ -358,7 +350,6 @@ public static class SimdExtensions var maxVec = minVec; int i = vectorSize; - // Process in vector chunks for (; i <= span.Length - vectorSize; i += vectorSize) { var vector = new Vector(span.Slice(i, vectorSize)); @@ -366,7 +357,6 @@ public static class SimdExtensions maxVec = Vector.Max(maxVec, vector); } - // Find min/max within vectors double min = minVec[0]; double max = maxVec[0]; for (int j = 1; j < vectorSize; j++) @@ -375,7 +365,6 @@ public static class SimdExtensions if (maxVec[j] > max) max = maxVec[j]; } - // Process remaining elements for (; i < span.Length; i++) { if (span[i] < min) min = span[i]; diff --git a/lib/core/tbar/tbar.cs b/lib/core/tbar/tbar.cs index 60bb4853..0175fedf 100644 --- a/lib/core/tbar/tbar.cs +++ b/lib/core/tbar/tbar.cs @@ -28,8 +28,8 @@ public readonly struct TBar : IEquatable // Computed properties (calculated on demand, no storage overhead) public double HL2 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (High + Low) * 0.5; } public double OC2 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (Open + Close) * 0.5; } - public double OHL3 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (Open + High + Low) / 3.0; } - public double HLC3 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (High + Low + Close) / 3.0; } + public double OHL3 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (Open + High + Low) * 0.333333333333333333; } + public double HLC3 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (High + Low + Close) * 0.333333333333333333; } public double OHLC4 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (Open + High + Low + Close) * 0.25; } public double HLCC4 { [MethodImpl(MethodImplOptions.AggressiveInlining)] get => (High + Low + Close + Close) * 0.25; } diff --git a/lib/core/tbarseries/tbarseries.cs b/lib/core/tbarseries/tbarseries.cs index 02fc022b..a54f5b02 100644 --- a/lib/core/tbarseries/tbarseries.cs +++ b/lib/core/tbarseries/tbarseries.cs @@ -11,7 +11,6 @@ namespace QuanTAlib; /// public class TBarSeries : IReadOnlyList { - // Internal storage: SoA layout protected readonly List _t = new(); protected readonly List _o = new(); protected readonly List _h = new(); @@ -22,7 +21,6 @@ public class TBarSeries : IReadOnlyList public string Name { get; set; } = "Bar"; public event Action? Pub; - // Public properties are Views into the main data public TSeries Open { get; } public TSeries High { get; } public TSeries Low { get; } @@ -38,7 +36,6 @@ public class TBarSeries : IReadOnlyList public TBarSeries() { - // Initialize views sharing the same Time list but different Value lists Open = new TSeries(_t, _o) { Name = "Open" }; High = new TSeries(_t, _h) { Name = "High" }; Low = new TSeries(_t, _l) { Name = "Low" }; @@ -46,9 +43,6 @@ public class TBarSeries : IReadOnlyList Volume = new TSeries(_t, _v) { Name = "Volume" }; } - /// - /// Constructor with capacity hint to avoid List growth overhead. - /// public TBarSeries(int capacity) { _t = new List(capacity); @@ -58,7 +52,6 @@ public class TBarSeries : IReadOnlyList _c = new List(capacity); _v = new List(capacity); - // Initialize views sharing the same Time list but different Value lists Open = new TSeries(_t, _o) { Name = "Open" }; High = new TSeries(_t, _h) { Name = "High" }; Low = new TSeries(_t, _l) { Name = "Low" }; diff --git a/lib/core/tseries/tseries.cs b/lib/core/tseries/tseries.cs index 9723f5ff..38dd1651 100644 --- a/lib/core/tseries/tseries.cs +++ b/lib/core/tseries/tseries.cs @@ -11,36 +11,25 @@ namespace QuanTAlib; /// public class TSeries : IReadOnlyList { - // Internal storage: SoA layout - // We use List for dynamic sizing but access internal arrays via CollectionsMarshal for speed protected readonly List _t; protected readonly List _v; public string Name { get; set; } = "Data"; - // Event optimization: Use Action to avoid EventArgs allocation - // Note: Events are generally discouraged in the hot path of this high-perf design, - // but kept for compatibility/chaining. public event Action? Pub; - public TSeries() + public TSeries() { _t = new List(); _v = new List(); } - /// - /// Constructor with capacity hint to avoid List growth overhead. - /// - public TSeries(int capacity) + public TSeries(int capacity) { _t = new List(capacity); _v = new List(capacity); } - /// - /// Constructor for wrapping existing lists (e.g. from TBarSeries). - /// public TSeries(List time, List values) { _t = time; @@ -105,7 +94,6 @@ public class TSeries : IReadOnlyList } else { - // Update last bar int lastIdx = _v.Count - 1; _t[lastIdx] = value.Time; _v[lastIdx] = value.Value; @@ -129,7 +117,7 @@ public class TSeries : IReadOnlyList foreach (var v in values) { Add(new TValue(t, v), isNew: true); - t += TimeSpan.TicksPerMinute; // Dummy time increment + t += TimeSpan.TicksPerMinute; } } diff --git a/lib/core/tvalue/tvalue.cs b/lib/core/tvalue/tvalue.cs index f39ee3e1..01bfca5a 100644 --- a/lib/core/tvalue/tvalue.cs +++ b/lib/core/tvalue/tvalue.cs @@ -9,19 +9,9 @@ namespace QuanTAlib; [SkipLocalsInit] public readonly struct TValue : IEquatable { - /// - /// Time in ticks (UTC). - /// public readonly long Time; - - /// - /// The value. - /// public readonly double Value; - /// - /// Convenience property to get DateTime from Ticks. - /// public DateTime AsDateTime => new(Time, DateTimeKind.Utc); [MethodImpl(MethodImplOptions.AggressiveInlining)] diff --git a/quantower/IndicatorExtensions.Tests.cs b/quantower/IndicatorExtensions.Tests.cs new file mode 100644 index 00000000..5970779d --- /dev/null +++ b/quantower/IndicatorExtensions.Tests.cs @@ -0,0 +1,148 @@ +using Xunit; +using TradingPlatform.BusinessLayer; +using System.Drawing; +using System.Reflection; + +namespace QuanTAlib.Tests; + +public class IndicatorExtensionsTests +{ + private class TestIndicator : Indicator + { + public TestIndicator() + { + Name = "Test Indicator"; + } + } + + private class TestCoordinatesConverter : ICoordinatesConverter + { + private readonly DateTime _time; + public TestCoordinatesConverter(DateTime time) => _time = time; + + public DateTime GetTime(int x) => _time; + public double GetChartX(DateTime time) => 0; + public double GetChartY(double value) => 0; + } + + [Fact] + public void DataSourceInputAttribute_HasCorrectDefaults() + { + var attr = new IndicatorExtensions.DataSourceInputAttribute(); + + Assert.Equal("Data source", attr.Name); + Assert.Equal(20, attr.SortIndex); + Assert.NotNull(attr.Variants); + Assert.NotEmpty(attr.Variants); + } + + [Fact] + public void GetInputValue_ReturnsCorrectValues_ForSourceTypes() + { + var indicator = new TestIndicator(); + var now = DateTime.UtcNow; + + // Open=100, High=110, Low=90, Close=105, Volume=1000 + indicator.HistoricalData.AddBar(now, 100, 110, 90, 105, 1000); + + // Ensure Count is updated (mock implementation detail) + // The mock HistoricalData.Count reflects added items. + // Indicator.Count => HistoricalData.Count. + + var args = new UpdateArgs(UpdateReason.NewBar); + + // Test each SourceType + Assert.Equal(100, IndicatorExtensions.GetInputValue(indicator, args, SourceType.Open).Value); + Assert.Equal(110, IndicatorExtensions.GetInputValue(indicator, args, SourceType.High).Value); + Assert.Equal(90, IndicatorExtensions.GetInputValue(indicator, args, SourceType.Low).Value); + Assert.Equal(105, IndicatorExtensions.GetInputValue(indicator, args, SourceType.Close).Value); + + // HL2 = (110 + 90) / 2 = 100 + Assert.Equal(100, IndicatorExtensions.GetInputValue(indicator, args, SourceType.HL2).Value); + + // OC2 = (100 + 105) / 2 = 102.5 + Assert.Equal(102.5, IndicatorExtensions.GetInputValue(indicator, args, SourceType.OC2).Value); + + // OHL3 = (100 + 110 + 90) / 3 = 100 + Assert.Equal(100, IndicatorExtensions.GetInputValue(indicator, args, SourceType.OHL3).Value); + + // HLC3 = (110 + 90 + 105) / 3 = 101.666... + Assert.Equal(101.66666666666667, IndicatorExtensions.GetInputValue(indicator, args, SourceType.HLC3).Value, 5); + + // OHLC4 = (100 + 110 + 90 + 105) / 4 = 101.25 + Assert.Equal(101.25, IndicatorExtensions.GetInputValue(indicator, args, SourceType.OHLC4).Value); + + // HLCC4 = (110 + 90 + 105 + 105) / 4 = 102.5 + Assert.Equal(102.5, IndicatorExtensions.GetInputValue(indicator, args, SourceType.HLCC4).Value); + } + + [Fact] + public void GetInputBar_ReturnsCorrectBar() + { + var indicator = new TestIndicator(); + var now = DateTime.UtcNow; + + indicator.HistoricalData.AddBar(now, 100, 110, 90, 105, 1000); + var args = new UpdateArgs(UpdateReason.NewBar); + + var bar = IndicatorExtensions.GetInputBar(indicator, args); + + Assert.Equal(now, bar.AsDateTime); + Assert.Equal(100, bar.Open); + Assert.Equal(110, bar.High); + Assert.Equal(90, bar.Low); + Assert.Equal(105, bar.Close); + Assert.Equal(1000, bar.Volume); + } + + [Fact] + public void PaintMethods_DoNotThrow_WithValidGraphics() + { + // This test attempts to verify that paint methods don't crash. + // It requires System.Drawing.Common to be functional. + + if (!System.Runtime.InteropServices.RuntimeInformation.IsOSPlatform(System.Runtime.InteropServices.OSPlatform.Windows)) + { + // Skip on non-Windows if System.Drawing is not fully supported (GDI+) + return; + } + + using var bitmap = new Bitmap(100, 100); + using var graphics = Graphics.FromImage(bitmap); + + var indicator = new TestIndicator(); + indicator.CurrentChart = new MockChart(); + + // Add some data + var now = new DateTime(2024, 1, 1, 12, 0, 0, DateTimeKind.Utc); + for (int i = 0; i < 20; i++) + { + indicator.HistoricalData.AddBar(now.AddMinutes(i), 100, 110, 90, 105); + } + + // Setup converter to return a time that exists in our data (e.g. the middle bar) + // We added bars at now, now+1min, ..., now+19min. + // Let's return now+10min. + var validTime = now.AddMinutes(10); + indicator.CurrentChart.MainWindow.CoordinatesConverter = new TestCoordinatesConverter(validTime); + + var args = new PaintChartEventArgs(graphics, new Rectangle(0, 0, 100, 100)); + using var pen = new Pen(Color.Red); + + // Test PaintHLine + IndicatorExtensions.PaintHLine(indicator, args, 100, pen); + + // Test PaintSmoothCurve + var series = new LineSeries("Test", Color.Blue, 1, LineStyle.Solid); + for (int i = 0; i < 20; i++) series.AddValue(); // Fill with NaNs or values + for (int i = 0; i < 20; i++) series.SetValue(100 + i, i); // Set some values + + IndicatorExtensions.PaintSmoothCurve(indicator, args, series, 0); + + // Test PaintHistogram + IndicatorExtensions.PaintHistogram(indicator, args, series, 0); + + // Test DrawText + IndicatorExtensions.DrawText(indicator, args, "Test Text"); + } +} diff --git a/quantower/TrimaIndicator.Tests.cs b/quantower/TrimaIndicator.Tests.cs new file mode 100644 index 00000000..a5763249 --- /dev/null +++ b/quantower/TrimaIndicator.Tests.cs @@ -0,0 +1,143 @@ +using Xunit; +using TradingPlatform.BusinessLayer; + +namespace QuanTAlib.Tests; + +public class TrimaIndicatorTests +{ + [Fact] + public void TrimaIndicator_Constructor_SetsDefaults() + { + var indicator = new TrimaIndicator(); + + Assert.Equal(10, indicator.Period); + Assert.Equal(SourceType.Close, indicator.Source); + Assert.True(indicator.ShowColdValues); + Assert.Equal("TRIMA - Triangular Moving Average", indicator.Name); + Assert.False(indicator.SeparateWindow); + Assert.True(indicator.OnBackGround); + } + + [Fact] + public void TrimaIndicator_MinHistoryDepths_EqualsPeriod() + { + var indicator = new TrimaIndicator { Period = 20 }; + + Assert.Equal(20, indicator.MinHistoryDepths); + Assert.Equal(20, ((IWatchlistIndicator)indicator).MinHistoryDepths); + } + + [Fact] + public void TrimaIndicator_ShortName_IncludesPeriodAndSource() + { + var indicator = new TrimaIndicator { Period = 15 }; + + Assert.Contains("TRIMA", indicator.ShortName); + Assert.Contains("15", indicator.ShortName); + } + + [Fact] + public void TrimaIndicator_SourceCodeLink_IsValid() + { + var indicator = new TrimaIndicator(); + + Assert.Contains("github.com", indicator.SourceCodeLink); + Assert.Contains("Trima.Quantower.cs", indicator.SourceCodeLink); + } + + [Fact] + public void TrimaIndicator_Initialize_CreatesInternalTrima() + { + var indicator = new TrimaIndicator { Period = 10 }; + + // Initialize should not throw + indicator.Initialize(); + + // After init, line series should exist + Assert.Single(indicator.LinesSeries); + } + + [Fact] + public void TrimaIndicator_ProcessUpdate_HistoricalBar_ComputesValue() + { + var indicator = new TrimaIndicator { Period = 3 }; + indicator.Initialize(); + + // Add historical data + var now = DateTime.UtcNow; + indicator.HistoricalData.AddBar(now, 100, 105, 95, 102); + + // Process update + var args = new UpdateArgs(UpdateReason.HistoricalBar); + indicator.ProcessUpdate(args); + + // Line series should have a value + Assert.Equal(1, indicator.LinesSeries[0].Count); + Assert.True(double.IsFinite(indicator.LinesSeries[0].GetValue(0))); + } + + [Fact] + public void TrimaIndicator_MultipleUpdates_ProducesCorrectTrimaSequence() + { + var indicator = new TrimaIndicator { Period = 3 }; + indicator.Initialize(); + + var now = DateTime.UtcNow; + double[] closes = { 100, 102, 104, 103, 105 }; + + foreach (var close in closes) + { + indicator.HistoricalData.AddBar(now, close, close + 2, close - 2, close); + indicator.ProcessUpdate(new UpdateArgs(UpdateReason.HistoricalBar)); + now = now.AddMinutes(1); + } + + // All values should be finite + for (int i = 0; i < closes.Length; i++) + { + Assert.True(double.IsFinite(indicator.LinesSeries[0].GetValue(closes.Length - 1 - i))); + } + + // TRIMA is smoothed, so check last value is reasonable + double lastTrima = indicator.LinesSeries[0].GetValue(0); + Assert.True(lastTrima >= 100 && lastTrima <= 106); + } + + [Fact] + public void TrimaIndicator_DifferentSourceTypes_Work() + { + var sources = new[] { SourceType.Open, SourceType.High, SourceType.Low, SourceType.Close, SourceType.HL2, SourceType.HLC3 }; + + foreach (var source in sources) + { + var indicator = new TrimaIndicator { Period = 3, Source = source }; + indicator.Initialize(); + + var now = DateTime.UtcNow; + indicator.HistoricalData.AddBar(now, 100, 110, 90, 105); + indicator.ProcessUpdate(new UpdateArgs(UpdateReason.HistoricalBar)); + + Assert.True(double.IsFinite(indicator.LinesSeries[0].GetValue(0)), + $"Source {source} should produce finite value"); + } + } + + [Fact] + public void TrimaIndicator_Period_CanBeChanged() + { + var indicator = new TrimaIndicator { Period = 5 }; + Assert.Equal(5, indicator.Period); + + indicator.Period = 20; + Assert.Equal(20, indicator.Period); + Assert.Equal(20, indicator.MinHistoryDepths); + } + + [Fact] + public void TrimaIndicator_DescriptionIsSet() + { + var indicator = new TrimaIndicator(); + + Assert.Contains("Triangular", indicator.Description); + } +}