using System.Runtime.CompilerServices;
namespace QuanTAlib;
///
/// BLMA: Blackman Moving Average
///
///
/// Window-based MA using Blackman coefficients (a0=0.42, a1=0.5, a2=0.08).
/// Minimizes spectral leakage with smooth taper to zero at edges.
///
/// Calculation: W_i = 0.42 - 0.5×cos(2πi/(n-1)) + 0.08×cos(4πi/(n-1)).
///
/// Detailed documentation
[SkipLocalsInit]
public sealed class Blma : AbstractBase
{
private readonly int _period;
private readonly RingBuffer _buffer;
private readonly double[] _weights;
private readonly double _weightSum;
private readonly TValuePublishedHandler _handler;
private ITValuePublisher? _source;
private int _disposed;
public override bool IsHot => _buffer.Count >= _period;
public Blma(int period)
{
if (period < 1)
{
throw new ArgumentOutOfRangeException(nameof(period), "Period must be greater than 0");
}
_period = period;
Name = $"Blma({period})";
WarmupPeriod = period;
_buffer = new RingBuffer(period);
_weights = new double[period];
_handler = Handle;
// Pre-calculate weights for the full period
_weightSum = CalculateWeights(period, _weights);
}
public Blma(ITValuePublisher source, int period) : this(period)
{
_source = source ?? throw new ArgumentNullException(nameof(source));
_source.Pub += _handler;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private void Handle(object? sender, in TValueEventArgs args)
{
Update(args.Value, args.IsNew);
}
///
/// Disposes the Blma instance, unsubscribing from the source publisher if subscribed.
/// This method is idempotent and thread-safe.
///
protected override void Dispose(bool disposing)
{
if (Interlocked.CompareExchange(ref _disposed, 1, 0) == 0 && _source != null)
{
_source.Pub -= _handler;
_source = null;
}
base.Dispose(disposing);
}
public override void Reset()
{
_buffer.Clear();
Last = default;
}
public override void Prime(ReadOnlySpan source, TimeSpan? step = null)
{
TimeSpan increment = step ?? TimeSpan.FromMilliseconds(1);
DateTime time = DateTime.UtcNow;
foreach (var value in source)
{
Update(new TValue(time, value));
time = time.Add(increment);
}
}
public void Prime(ReadOnlySpan source)
{
foreach (var value in source)
{
Update(value);
}
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public override TValue Update(TValue input, bool isNew = true)
{
// Handle NaN/Infinity - return last result without changing state
double val = input.Value;
if (!double.IsFinite(val))
{
return Last;
}
_buffer.Add(val, isNew);
double result;
if (_buffer.Count < _period)
{
// During warmup, calculate weights dynamically for the current count
int count = _buffer.Count;
if (count == 1)
{
result = val;
}
else
{
Span currentWeights = stackalloc double[count];
double currentWeightSum = CalculateWeights(count, currentWeights);
result = ComputeWeightedAverage(
currentWeightSum,
CalculateWeightedSum(_buffer, currentWeights),
_buffer.Average);
}
}
else
{
// Full period, use pre-calculated weights
result = ComputeWeightedAverage(
_weightSum,
CalculateWeightedSum(_buffer, _weights),
_buffer.Average);
}
var tValue = new TValue(input.Time, result);
Last = tValue;
PubEvent(tValue, isNew);
return tValue;
}
public override TSeries Update(TSeries source)
{
var result = new TSeries();
int len = source.Count;
// Use ArrayPool for large allocations
double[]? rented = len > 256 ? System.Buffers.ArrayPool.Shared.Rent(len) : null;
scoped Span output = rented != null ? rented.AsSpan(0, len) : stackalloc double[len];
try
{
BatchCore(source.Values, output, _period);
for (int i = 0; i < len; i++)
{
result.Add(new TValue(source[i].Time, output[i]));
}
}
finally
{
if (rented != null)
{
System.Buffers.ArrayPool.Shared.Return(rented);
}
}
// Restore state by replaying last Period bars
// This ensures the indicator is ready for subsequent streaming updates
Reset();
int start = Math.Max(0, len - _period);
for (int i = start; i < len; i++)
{
Update(source[i]);
}
return result;
}
///
/// Computes weighted average with fallback for zero weight sum.
///
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static double ComputeWeightedAverage(double weightSum, double weightedSum, double fallbackAverage)
{
return Math.Abs(weightSum) < double.Epsilon ? fallbackAverage : weightedSum / weightSum;
}
private static double CalculateWeights(int n, Span weights)
{
if (n == 1)
{
weights[0] = 1.0;
return 1.0;
}
double totalWeight = 0;
double invNMinus1 = 1.0 / (n - 1);
const double pi2 = 2.0 * Math.PI;
const double pi4 = 4.0 * Math.PI;
// Blackman window coefficients
const double a0 = 0.42;
const double a1 = 0.5;
const double a2 = 0.08;
for (int i = 0; i < n; i++)
{
double ratio = i * invNMinus1;
// Use FMA: a0 - a1*cos1 + a2*cos2 = a0 + FMA(-a1, cos1, a2*cos2)
double cos1 = Math.Cos(pi2 * ratio);
double cos2 = Math.Cos(pi4 * ratio);
double w = a0 + Math.FusedMultiplyAdd(-a1, cos1, a2 * cos2);
weights[i] = w;
totalWeight += w;
}
return totalWeight;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static double CalculateWeightedSum(RingBuffer buffer, ReadOnlySpan weights)
{
int start = buffer.StartIndex;
int count = buffer.Count;
int capacity = buffer.Capacity;
if (start + count <= capacity)
{
return buffer.InternalBuffer.Slice(start, count).DotProduct(weights);
}
int firstPartLength = capacity - start;
int secondPartLength = count - firstPartLength;
double sum1 = buffer.InternalBuffer.Slice(start, firstPartLength).DotProduct(weights[..firstPartLength]);
double sum2 = buffer.InternalBuffer.Slice(0, secondPartLength).DotProduct(weights[firstPartLength..]);
return sum1 + sum2;
}
///
/// Calculates BLMA values for a TSeries and returns both results and a primed indicator.
///
///
/// Core batch calculation of BLMA values using spans (high-performance batch API).
///
private static void BatchCore(ReadOnlySpan source, Span destination, int period)
{
if (period < 1)
{
throw new ArgumentOutOfRangeException(nameof(period), "Period must be greater than 0");
}
if (destination.Length < source.Length)
{
throw new ArgumentOutOfRangeException(nameof(destination), $"Destination length must be at least {source.Length}.");
}
// Pre-calculate weights for full period
Span weights = period <= 256 ? stackalloc double[period] : new double[period];
double weightSum = CalculateWeights(period, weights);
// Buffer for warmup weights to avoid stackalloc in loop
Span warmupWeightsBuffer = period <= 256 ? stackalloc double[period] : new double[period];
// Handle NaN via last-valid-value substitution
double lastValid = double.NaN;
for (int i = 0; i < source.Length; i++)
{
if (double.IsFinite(source[i]))
{
lastValid = source[i];
break;
}
}
for (int i = 0; i < source.Length; i++)
{
double val = source[i];
if (!double.IsFinite(val))
{
val = double.IsNaN(lastValid) ? 0 : lastValid;
}
else
{
lastValid = val;
}
int count = Math.Min(i + 1, period);
if (count < period)
{
// Warmup: dynamic weights
if (count == 1)
{
destination[i] = val;
}
else
{
Span currentWeights = warmupWeightsBuffer.Slice(0, count);
double currentWeightSum = CalculateWeights(count, currentWeights);
double sum = 0;
for (int j = 0; j < count; j++)
{
int srcIdx = i - count + 1 + j;
double srcVal = source[srcIdx];
if (!double.IsFinite(srcVal))
{
srcVal = lastValid;
}
sum += srcVal * currentWeights[j];
}
double avg = 0;
for (int j = 0; j < count; j++)
{
int srcIdx = i - count + 1 + j;
double srcVal = source[srcIdx];
if (!double.IsFinite(srcVal))
{
srcVal = lastValid;
}
avg += srcVal;
}
avg /= count;
destination[i] = ComputeWeightedAverage(currentWeightSum, sum, avg);
}
}
else
{
// Full period
double sum = 0;
double avg = 0;
for (int j = 0; j < period; j++)
{
int srcIdx = i - period + 1 + j;
double srcVal = source[srcIdx];
if (!double.IsFinite(srcVal))
{
srcVal = lastValid;
}
sum += srcVal * weights[j];
avg += srcVal;
}
avg /= period;
destination[i] = ComputeWeightedAverage(weightSum, sum, avg);
}
}
}
///
/// Batch calculates BLMA values for a TSeries.
///
public static TSeries Batch(TSeries source, int period)
{
var indicator = new Blma(period);
return indicator.Update(source);
}
///
/// Batch calculates BLMA values using spans.
///
public static void Batch(ReadOnlySpan source, Span destination, int period)
{
BatchCore(source, destination, period);
}
public static (TSeries Results, Blma Indicator) Calculate(TSeries source, int period)
{
var indicator = new Blma(period);
var results = indicator.Update(source);
return (results, indicator);
}
}