Extrapolator - sincMA

//+------------------------------------------------------------------+
//|                                        Extrapolator - sincMA.mq4 |
//|                                          Copyright © 2009, gpwr. |
//+------------------------------------------------------------------+
#property copyright "Copyright © 2009, gpwr."
#property indicator_chart_window
#property indicator_buffers 2
#property indicator_color1 Blue
#property indicator_color2 Red
#property indicator_width1 2
#property indicator_width2 2

//Global constants
#define pi 3.141592653589793238462643383279502884197169399375105820974944592

//Input parameters
extern int    LastBar  =0;   //Last bar in the past data
extern int    PastBars =512; //Number of past bars. 0: all data.
extern int    FutBars  =30;  //Number of predicted future bars. Must be less than PastBars
extern int    trun     =3;   //Number of the last wavelets to truncate.
extern double Q        =1.0; //Damping factor. Higher Q gives more ripple. Q=1 gives sinc function. 

//Global variables
int n,nf,lb,pwrMax; 

//Indicator buffers
double pv[],fv[];

//Initialize---------------------------------------------------------+
int init()
  {
   // Initialize number of fitted values and the max number of sinc pulses
   lb=LastBar;
   nf=FutBars;
   if(PastBars<=0) PastBars=Bars-lb;
   n=1;
   int pwr=0;
   while(n<=PastBars)
   {
      pwr++;
      n*=2;
   }
   pwr--;
   n/=2;
   if(nf>n) nf=n;
   if(trun<=0) trun=0;
   if(trun>=pwr) trun=pwr;
   pwrMax=pwr-trun;
   
   //Initialize indicators
   IndicatorBuffers(2);
   SetIndexBuffer(0,pv);
   SetIndexBuffer(1,fv);
   SetIndexStyle(0,DRAW_LINE,STYLE_SOLID,2);
   SetIndexStyle(1,DRAW_LINE,STYLE_SOLID,2);
   SetIndexShift(0,-lb);//past data vector 0..n; 0 corresponds to bar=lb
   SetIndexShift(1,nf-lb);//future data vector i=0..nf; nf corresponds to bar=lb
   IndicatorShortName("Extrapolator - sincMA");  
   return(0);
  }

int deinit(){return(0);} 

//START--------------------------------------------------------------+
int start()
{
   // Initialize the indicator buffers with straight lines
   for(int i=0;i<=n;i++)
   {
      pv[i]=Close[lb]+(Close[n+lb]-Close[lb])*i/n;
      if(i<=nf) fv[i]=Close[lb]-(Close[n+lb]-Close[lb])*i/n;
   }
   
   // Sample the data and fit sinc pulses to it
   double dp[],df[];
   ArrayResize(dp,n+1);
   ArrayResize(df,nf+1);
   int ss=2*n; // sampling step
   for(int p=1;p<=pwrMax;p++)
   {
      ss=ss/2;
      for(i=0;i<=n;i++)
      {
         dp[i]=0.0;
         if(i<=nf) df[i]=0.0;
         for(int j=1;j<=n/ss;j++)
         {
            double dif=Close[j*ss-ss/2+lb]-pv[j*ss-ss/2];
            dp[i]+=dif*wave(2*j-1-2.0*i/ss,Q);
            if(i<=nf) df[i]+=dif*wave(2*j-1+2.0*i/ss,Q);
         }
      }
      for(i=0;i<=n;i++)
      {
         pv[i]+=dp[i];
         if(i<=nf) fv[i]+=df[i];
      }
   }
   
   // Reorder fv[]
   for(i=0;i<nf/2;i++)
   {
      double tmp=fv[i];
      fv[i]=fv[nf-i];
      fv[nf-i]=tmp;
   }
   
   return(0);
}

//sinc function------------------------------------------------------+
double wave(double x, double Q)
{
   if(x==0.0) return(1.0);
   else return(MathSin(pi*x)/(pi*x)/(x*x+1.0)*(Q*x*x+1.0));
}
//+------------------------------------------------------------------+