#import <objc/Object.h>
#import "ART_STM.h"
#include "ART_Func.h"

#define ART_ACT_OFF -1

@implementation ART_STM
/*
float   *P,*U,*W,         // left Side
        *Q,*V,*X,         // right Side
        *old_W,*old_P,    // last Iteration state
        *I,*E;            // input at W and input at P

float   a,b,              // field constances
        d,                // learnrate
        theta;            // parameter of the non-linear function

int     i_min,i_max,      // iteration boundings
        act_node;         // the active node of a competionen fields
float   Act_Stable;       // activity stable criterion
*/


+ new
{
     self=[super new];
     return self;
}

- save
{
     return self;
}

- load
{
     return self;
}


- init:(int)Vec_dim
{
     int Heapspace=Vec_dim*sizeof(float);
     [super init];
     M=Vec_dim;
     zone=[self zone];
     P=NXZoneMalloc(zone,Heapspace);
     U=NXZoneMalloc(zone,Heapspace);
     W=NXZoneMalloc(zone,Heapspace);
     Q=NXZoneMalloc(zone,Heapspace);
     V=NXZoneMalloc(zone,Heapspace);
     X=NXZoneMalloc(zone,Heapspace);
     old_W=NXZoneMalloc(zone,Heapspace);
     old_P=NXZoneMalloc(zone,Heapspace);
     [self set_zero];
     [self set_a:10 b:10 theta:(1.1/(float)sqrt(Vec_dim)) d:1.0];
     [self set_NoLinFunc:ART_STM_SIGMO];
     act_node=ART_ACT_OFF;
     return self;
}
- free
{
     NXZoneFree(zone,P);
     NXZoneFree(zone,U);
     NXZoneFree(zone,W);
     NXZoneFree(zone,Q);
     NXZoneFree(zone,V);
     NXZoneFree(zone,X);
     NXZoneFree(zone,old_W);
     NXZoneFree(zone,old_P);     
     return [super free];
}

- set_a:(float)p1 b:(float)p2 theta:(float)p3 d:(float)p4
{
     a=p1;
     b=p2;
     theta=p3;
     d=p4;
     return self;
}
- set_NoLinFunc:(int)f
{
     noLinFunc=f;
     return self;
}

- set_zero
{
     int i;
     for (i=0;i<M;i++)
       P[i]=U[i]=W[i]=0.0;
     act_node=ART_ACT_OFF;
     return self;
}

- set_I:(float *)I_PTR and_E:(float *)E_PTR
{
     I=I_PTR;
     E=E_PTR;
     return self;
}

- setIterationmin:(int)min max:(int)max stable:(float)stable
{
     i_min=min;
     i_max=max;
     Act_stable=stable;
     return self;
}

- (float *)get_P { return P; }
- (float *)get_U { return U; }
- (float *)get_W { return W; }
- (float *)get_Q { return Q; }
- (float *)get_V { return V; }
- (float *)get_X { return X; }

- (float *)get_I { return I; }
- (float *)get_E { return E; }

- set_P:(float *)Vec
{
      P=Vec;
      return self;
}


// Special F2 competition method retrun the index of maximum of U
- (int)MaxAct_with:(int)nodes;
{
     int i;
     float vec_max=0.0;
     act_node=ART_ACT_OFF;
     for (i=0; i<nodes; i++)
       vec_max = (vec_max < U[i] ? U[i] : vec_max);
    if (vec_max == 0.0) return (act_node);
     for (i=0; i<nodes && U[i]<vec_max ; i++);
     act_node=i;
     return (act_node) ;

}
-(int)activeNode
{
     return(act_node);
}


// STM Object Iteration and Performance Method
- (int)iterate
{
     int i,cntr=0;
     float ptmp,wtmp,vtmp;
     float v_d = 1.0;

     do {
       // Pass I of the iteration
       for (i=0; i<M; i++)  {
         old_W[i]=W[i];
         old_P[i]=P[i];
         W[i]=I[i]+a*U[i];
         P[i]=U[i]+d*E[i];
         }

       ptmp=L2_norm(P,M);
       wtmp=L2_norm(W,M);
       ptmp = (ptmp > 0.1) ? ptmp : 1.0;
       wtmp = (wtmp > 0.1) ? wtmp : 1.0;

       // Pass II of the Iteration
       switch(noLinFunc) {
         case ART_STM_SIGMO:  
           for (i=0; i<M; i++)  {
             Q[i]=P[i]/ptmp;
             X[i]=W[i]/wtmp;
             V[i]=sigmo(X[i],theta,0.001)+b*sigmo(Q[i],theta,0.001);
           }
           break;
         case ART_STM_GROSS:  
           for (i=0; i<M; i++)  {
             Q[i]=P[i]/ptmp;
             X[i]=W[i]/wtmp;
             V[i]=fgross(X[i],theta)+b*fgross(Q[i],theta);
           }
           break;
         }  



       vtmp=L2_norm(V,M);
       vtmp = (vtmp > 0.1) ? vtmp : 1.0;

       // Pass III of the iteration
       for (i=0; i<M; i++)
         U[i]=V[i]/vtmp;

       // Check for Stability by means of the absolute Vectordiff.
       v_d=vec_diff(old_W,W,M)+vec_diff(old_P,P,M);

       // Incrment the iterationcounter
       cntr++;
       }
     // Do Iterate again
     while (
       // 1. W and P aren't stable and iteraion didn't excced i_max
       ( v_d > Act_stable && cntr < i_max )
       // 2. iteration are less then i_min
       || cntr < i_min
     );
     return (int) cntr;
}

@end





