view release on metacpan or  search on metacpan

lib/App/Chart/Series/Derived/EMAx2.pm  view on Meta::CPAN

use Locale::TextDomain 1.17; # for __p()
use Locale::TextDomain ('App-Chart');

use base 'App::Chart::Series::Indicator';
use App::Chart::Series::Derived::EMA;

sub longname   { __('EMA of EMA') }
sub shortname  { __('EMAofEMA') }
sub manual     { __p('manual-node','EMA of EMA') }

use constant
  { type       => 'average',
    priority   => -10,
    parameter_info => [ { name     => __('Days'),
                          key      => 'ema2_days',
                          type     => 'float',
                          minimum  => 1,
                          default  => 20,
                          decimals => 0,
                          step     => 1 } ],
  };

sub new {
  my ($class, $parent, $N) = @_;

  $N //= parameter_info()->[0]->{'default'};
  ($N > 0) or croak "EMA2 bad N: $N";

  return $class->SUPER::new
    (parent     => $parent,
     parameters => [ $N ],
     N          => $N,
     arrays     => { values => [] },
     array_aliases => { });
}

# Return a procedure which calculates an EMA of EMA with given $N period
# smoothing (on both).
#
# Each call $proc->($value) enters a new value into the window, and the
# return is the EMAofEMA up to (and including) that value.
#
# An EMA of EMA is in theory influenced by all preceding data, but
# warmup_count() below is designed to determine a warmup count.
#
sub proc {
  my ($class_or_self, $N) = @_;
  my $ema_proc = App::Chart::Series::Derived::EMA->proc ($N);
  my $ema2_proc = App::Chart::Series::Derived::EMA->proc ($N);
  return sub { $ema2_proc->($ema_proc->($_[0])) };
}

# By priming an EMA-2 pro with warmup_count() many values, the next call
# will have an omitted weight of no more than 0.1% of the total.  Omitting
# 0.1% should be negligable, unless past values are ridiculously bigger than
# recent ones.
#
# The implementation here does a binary search for the first k satisfying
# R(k)<=0.001, so it's only moderately fast.  Perhaps there'd be a direct
# closed-form solution to the equation R(k)=0.001 below.  The inverse of
# k*f^k is close to the Lambert W-function.  Is there an easy formula for
# that?
#
sub warmup_count {
  my ($self_or_class, $N) = @_;

  if ($N <= 1) { return 0; }
  my $f = App::Chart::Series::Derived::EMA::N_to_f ($N);
  return bsearch_first_true
    (sub {
       my ($i) = @_;
       return (ema2_omitted($f,$i)
               <= App::Chart::Series::Derived::EMA::WARMUP_OMITTED_FRACTION);
     },
     $N);
}

# ema2_omitted() returns the fraction (between 0 and 1) of weight omitted
# by stopping an EMA of EMA at the f^k term, which means the first k+1
# terms.
#
# The total weight up to that term is,
#
#    W(k) = (1-f)^2 * ( 1 + 2f + 3f^2 + 4f^3 + ... + (k+1)*f^k )
#
# Multiplying (1-f)^2 through leads to the middle terms cancelling, leaving
# just the 1 at the start, and two terms at the end,
#
#    W(k) = 1 + (-2*(k+1) + k) * f^(k+1)
#             + (k+1)          * f^(k+2)
#
# The omitted part is 1 - W(k), so the "1 +" is dropped and the rest
# negated.  The terms simplify to
#
#    R(k) = f^(k+1) * (k+2 - f * (k+1))
#
# See devel/ema-omitted.pl for automated checking of this calculation.
#
sub ema2_omitted {
  my ($f, $k) = @_;
  return $f**($k+1) * ($k+2 - $f*($k+1));
}

sub bsearch_first_true {
  my ($pred, $incr) = @_;
  $incr = POSIX::ceil ($incr);

  my $lo = 0;
  my $hi = $incr;

  # search upwards
  until ($pred->($hi)) {
    $lo = $hi + 1;
    $hi *= 2;
  }

  # at this point $pred->($lo) is unknown, $pred->($hi) is true
  if ($pred->($lo)) { return $lo; }

  # binary search loop, with $pred->($lo) false, $pred->($hi) true
  for (;;) {