/******************************************************************************/
/*                                                                            */
/*                      Fuzzy set                                             */
/*             2001.08.03 by CHEN Ching-Han                                   */
/*                                                                            */
/******************************************************************************/

#if !defined(FZYSET_H)
#define FZYSET_H
#include <fstream.h>
#include <math.h>
#include <string>
#include <iomanip.h>
#include "array.h"
#include "fzylv.h"
class FuzzySet
{
public:
   FuzzySet (){ nb = 0 ; }
   void Initialize ( int npts , float *xpts , float *ypts ) ;
   void Initialize ( f1D &xpts , f1D &ypts ) ;
   void Initialize ( float x0,float y0,float x1,float y1 ,float x2,float y2) ;
   void Initialize ( float x0 , float y0 , float x1 , float y1 ,
              float x2 , float y2 , float x3 , float y3 ) ;
   void Initialize ( const FuzzySet& s ) ;
   FuzzySet& operator= ( const FuzzySet& s ) ;
   ~FuzzySet () ;

// defuzzification
   float centroid () ;

// conjunction and disjunction operation for fuzzy inferencing
// if conj=1 =>conjunction, else disjunction
   void AND(const FuzzySet& s){cj(s,1);}
   void OR(const FuzzySet& s){cj(s,0);}
   void negate () ;
   void scale ( float fac ) ;

// utilities for fuzzy set
   int invalid(){return !nb ;}
   float membership ( float pt ) ;

   void Print(ostream &out);
   void SetName(char *s){name = s;}
   void Labeling(int ndx){label = ndx;}
   void Copy(LinguisticVariable &LV, int ndx_fzyset);
   void Copy(FuzzySet &fzy);
   char *name;
   int label;
protected:
   void cj ( const FuzzySet& s , int conj ) ;
   int nb ;    // Number of points in membership function def
   f1D x ;     // Domain points
   f1D y ;     // Corresponding function values
private:
void cj_out ( float x , float y , f1D &xn , f1D &yn , int &nn ) ;
void intsec ( float x1 , float y1 , float x2 , float y2 ,
              float y3 , float y4 , float &xint , float &yint ) ;
};

void FuzzySet::Initialize ( int npts , float *xpts , float *ypts )
{
   if (npts < 1){ nb = 0 ;  return ; }
   nb = npts;
   x.m = xpts ;
   y.m = ypts ;
}

void FuzzySet::Initialize ( f1D &xpts , f1D &ypts )
{
   if (xpts.nb != ypts.nb){ nb = 0 ;  return ; }
   nb = xpts.nb;
   x.m = xpts.m ;
   y.m = ypts.m ;
}

void FuzzySet::Initialize ( float x0 , float y0 , float x1 , float y1 ,
                     float x2 , float y2 , float x3 , float y3 )
{
   float xx[4], yy[4] ;
   xx[0] = x0 ;
   yy[0] = y0 ;
   nb = 1 ;
   if ((x1 != x0)  ||  (y1 != y0))
   {
     xx[nb] = x1 ;
     yy[nb] = y1 ;
     ++nb ;
   }
   if ((x2 != x1)  ||  (y2 != y1))
   {
      xx[nb] = x2 ;
      yy[nb] = y2 ;
      ++nb ;
   }
   if ((x3 != x2)  ||  (y3 != y2))
   {
      xx[nb] = x3 ;
      yy[nb] = y3 ;
      ++nb ;
   }
   x.Initialize(nb) ;
   y.Initialize(nb) ;
   for(int i=0;i<nb;i++)
   {
      x[i] = xx[i];
      y[i] = yy[i];
   }
}
void FuzzySet::Initialize(float x0,float y0,float x1,float y1,float x2,float y2)
{
   float xx[3], yy[3] ;
   xx[0] = x0 ;
   yy[0] = y0 ;
   nb = 1 ;
   if ((x1 != x0)  ||  (y1 != y0))
   {
     xx[nb] = x1 ;
     yy[nb] = y1 ;
     ++nb ;
   }
   if ((x2 != x1)  ||  (y2 != y1))
   {
      xx[nb] = x2 ;
      yy[nb] = y2 ;
      ++nb ;
   }
   x.Initialize(nb) ;
   y.Initialize(nb) ;
   for(int i=0;i<nb;i++)
   {
      x[i] = xx[i];
      y[i] = yy[i];
   }
}

void FuzzySet::Initialize ( const FuzzySet& s )
{
   if (nb != 0)
   {
     x.Finish();
     y.Finish();
   }
   nb = s.nb;
   x.Initialize(nb) ;
   y.Initialize(nb) ;
   for(int i=0;i<nb;i++)
   {
      x.m[i] = s.x.m[i];
      y.m[i] = s.y.m[i];
   }
   label = s.label;
   name = s.name;
}


FuzzySet::~FuzzySet ()
{
   if(nb!=0)
   {
      x.Finish();
      y.Finish();
   }
}

FuzzySet& FuzzySet::operator= ( const FuzzySet& s )
{
   if (s.nb==0)return *this;
   if (nb != 0)
   {
     x.Finish();
     y.Finish();
   }
   nb = s.nb;
   x.Initialize(nb) ;
   y.Initialize(nb) ;
   for(int i=0;i<nb;i++)
   {
      x.m[i] = s.x.m[i];
      y.m[i] = s.y.m[i];
   }
   label = s.label;
   name = s.name;
   return *this;
}

void FuzzySet::Print(ostream &out)
{
  out<<"Fuzzy Set "<<label<<endl;
  for(int i=0;i<nb;i++)out<<setprecision(2)<<x[i]<<"\t"<<setprecision(2)<<y[i]<<endl;
}

/*
--------------------------------------------------------------------------------

   Membership evaluator, invalid indicator, printer

--------------------------------------------------------------------------------
*/

float FuzzySet::membership ( float pt )
{
   if (nb==0) return 0.0 ;
   int lo, mid, hi ;
   float yy ;
   if (pt <= x[0]) return y[0] ;
   if (pt >= x[nb-1]) return y[nb-1] ;
   lo = 0 ;    // We will keep x[lo] strictly less than pt
   hi = nb-1 ;  // and x[hi] greater or equal to pt

   for (;;)
   {                           // Cuts interval in half each time
      mid = (lo + hi) / 2 ;    // Center of interval
      if (mid == lo)break ;    // Happens when lo and hi adjacent
      if (x[mid]<pt)lo=mid;    // Replace appropriate interval end with mid
      else  hi = mid ;
   }
   yy = (pt - x[hi-1]) / (x[hi] - x[hi-1]) * (y[hi] - y[hi-1]) ;
   return yy + y[hi-1] ;
}

/*
--------------------------------------------------------------------------------

   Negation, Inference by scaling

--------------------------------------------------------------------------------
*/

void FuzzySet::negate ()
{
   int i = nb ;
   while (i--) y[i] = 1.0 - y[i] ;
}

void FuzzySet::scale ( float fac )
{
   int i = nb ;
   while (i--) y[i] = fac * y[i] ;
}

/*
--------------------------------------------------------------------------------

   Conjunction, disjunction

--------------------------------------------------------------------------------
*/



void FuzzySet::cj ( const FuzzySet& s, int conj)
{
   if ((!nb)||(!s.nb)) return ;
   int i, nn, vertex0, vertex1, use_func_0 ;
   f1D xn, yn;
   float xl, xr, y0l, y0r, y1l, y1r, xint, yint ;
   float next_x0, next_y0, next_x1, next_y1, rightmost_x, frac ;

   xn.Initialize(2 * (nb + s.nb) ) ;
   yn.Initialize(2 * (nb + s.nb) ) ;
   conj = (conj != 0) ;
/*
   In the following, the suffix '0' refers to the current set,
   being modified.  '1' means the external, modifying set.
   The domain endpoints are always xl and xr.  The corresponding
   function values are y0l, y0r, y1l, y1r.

   Initialize the leftmost point.  If one of the domains starts
   lower than the other, imply a continued function value for the
   other so they can start together.
*/
   if (x.m[0] < s.x.m[0])
   {
      xl = x.m[0] ;
      vertex0 = 1 ;
      vertex1 = 0 ;
   }
   else if (x.m[0] > s.x.m[0])
   {
      xl = s.x.m[0] ;
      vertex0 = 0 ;
      vertex1 = 1 ;
   }
   else
   {
      xl = x.m[0] ;
      vertex0 = vertex1 = 1 ;
   }
   y0l = y.m[0] ;
   y1l = s.y.m[0] ;

   nn = 0 ;
   cj_out ( xl , (conj ^ (y0l>y1l)) ? y0l : y1l , xn , yn , nn ) ;

//   Also initialize the rightmost point so we end together.

   if (x.m[nb-1] >= s.x[s.nb-1]) rightmost_x = x.m[nb-1] ;
   else rightmost_x = s.x.m[s.nb-1] ;

//   At any time from now on, next_?? will be the next point in
//   each function.  We will advance by choosing the closest.

   if (vertex0 < nb)
   {
      next_x0 = x.m[vertex0] ;
      next_y0 = y.m[vertex0] ;
   }
   else
   {
      next_x0 = rightmost_x ;
      next_y0 = y.m[vertex0-1] ;
   }

   if (vertex1 < s.nb)
   {
      next_x1 = s.x.m[vertex1] ;
      next_y1 = s.y.m[vertex1] ;
   }
   else
   {
      next_x1 = rightmost_x ;
      next_y1 = s.y.m[vertex1-1] ;
   }

//   Main loop is here

   while ((vertex0 < nb) || (vertex1 < s.nb))
   {
      if (next_x0 < next_x1) use_func_0 = 1 ;
      else if (next_x1 < next_x0) use_func_0 = 0 ;
      else use_func_0 = (vertex0 < nb) ;
      if (use_func_0)
      {
         xr = next_x0 ;
         y0r = next_y0 ;
         if (next_x1 == xl) frac = 0.0 ;
         else frac = (xr - xl) / (next_x1 - xl) ;
         y1r = y1l + frac * (next_y1 - y1l) ;

         if (++vertex0 < nb) { next_x0 = x.m[vertex0] ; next_y0 = y.m[vertex0] ; }
         else next_x0 = rightmost_x ;
      }
      else
      {
         xr = next_x1 ;
         y1r = next_y1 ;
         if (next_x0 == xl) frac = 0.0 ;
         else frac = (xr - xl) / (next_x0 - xl) ;
         y0r = y0l + frac * (next_y0 - y0l) ;
         if (++vertex1 < s.nb)
         {
           next_x1 = s.x.m[vertex1];
           next_y1 = s.y.m[vertex1];
         }
         else next_x1 = rightmost_x ;
      }
//   If the functions cross, generate new point at intersection
      if ((xr > xl)  &&  ((y0l-y1l) * (y0r-y1r) < 0.0))
      {
         intsec ( xl , y0l , xr , y0r , y1l , y1r , xint , yint ) ;
         cj_out ( xint , yint , xn , yn , nn ) ;
      }
      cj_out ( xr , (conj ^ (y0r>y1r)) ? y0r : y1r , xn , yn , nn ) ;
      xl = xr ;
      y0l = y0r ;
      y1l = y1r ;
   } // Main loop

   x.Finish();
   y.Finish();
   nb=nn;
   x.Initialize(nb) ;
   y.Initialize(nb) ;
   for (i=0 ; i<nb ; i++)
   {
   	x.m[i]=xn.m[i];
    y.m[i]=yn.m[i];
   }
}

void FuzzySet::cj_out ( // save a resultant point
   float x ,         // x value of new point
   float y ,         // and y
   f1D &xn ,         // Array of x's
   f1D &yn ,         // and y's
   int &nn)           // Number in above
{
   float d ;
   if((nn>0)&&(fabs(x-xn.m[nn-1])<1.e-20)&&(fabs(y-yn.m[nn-1])<1.e-20))return;
   if (nn > 1)
   {
      d = xn.m[nn-1] * (y - yn.m[nn-2]) +
          yn.m[nn-1] * (xn.m[nn-2] - x) +
          x * yn.m[nn-2] - y * xn.m[nn-2] ;
      if (fabs(d) < 1.e-10) --nn ;
   }
   xn.m[nn] = x ;
   yn.m[nn] = y ;
   ++nn ;
}

//
// intersect two line segments
//

void FuzzySet::intsec ( float x1 , float y1 , float x2 , float y2 ,
                     float y3 , float y4 , float &xint , float &yint )
{
   float den = y1 - y2 - y3 + y4 ;
   xint = (x2 * y1 - x1 * y2 - x2 * y3 + x1 * y4) / den ;
   yint = (y1 * y4 - y2 * y3) / den ;
}

/*
--------------------------------------------------------------------------------

   centroid - Compute the centroid of a fuzzy set

--------------------------------------------------------------------------------
*/

float FuzzySet::centroid ()
{
   float t, numer, denom ;
   if (nb < 2) return 0.0 ;
   numer = denom = 0.0 ;
   for (int i=1 ; i<nb ; i++)
   {
      t = y.m[i-1] * (2.0 * x.m[i-1] + x.m[i]) + y.m[i] * (x.m[i-1] + 2.0 * x.m[i]) ;
      numer += t * (x.m[i] - x.m[i-1]) ;
      denom += (x.m[i] - x.m[i-1]) * (y.m[i-1] + y.m[i]) ;
   }
   return numer / (3.0 * denom) ;
}

void FuzzySet::Copy(LinguisticVariable &LV, int ndx_fzyset)
{
  Initialize(LV.FzySet[ndx_fzyset]);
}

void FuzzySet::Copy(FuzzySet &fzy)
{
  Initialize(fzy);
}
#endif