/*
//
//                  INTEL CORPORATION PROPRIETARY INFORMATION
//     This software is supplied under the terms of a license agreement or
//     nondisclosure agreement with Intel Corporation and may not be copied
//     or disclosed except in accordance with the terms of that agreement.
//          Copyright(c) 1999-2006 Intel Corporation. All Rights Reserved.
//
//     Intel(R) Integrated Performance Primitives Gaussian Mixture Sample for Windows*
//
//  By downloading and installing this sample, you hereby agree that the
//  accompanying Materials are being provided to you under the terms and
//  conditions of the End User License Agreement for the Intel(R) Integrated
//  Performance Primitives product previously accepted by you. Please refer
//  to the file ippEULA.rtf located in the root directory of your Intel(R) IPP
//  product installation for more information.
//
//
//     Gaussian Selection Class
//
*/


#include "select_gauss.h"
#include "prob_calc.h"

int CountArray[]={0,1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,1,2,2,3,2,3,3,4,2,3,3,4,3,4,
 4,5,1,2,2,3,2,3,3,4,2,3,3,4,3,4,4,5,2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,1,2,2,3,
 2,3,3,4,2,3,3,4,3,4,4,5,2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,2,3,3,4,3,4,4,5,3,4,
 4,5,4,5,5,6,3,4,4,5,4,5,5,6,4,5,5,6,5,6,6,7,1,2,2,3,2,3,3,4,2,3,3,4,3,4,4,5,
 2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,3,4,4,5,4,5,
 5,6,4,5,5,6,5,6,6,7,2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,3,4,4,5,4,5,5,6,4,5,5,6,
 5,6,6,7,3,4,4,5,4,5,5,6,4,5,5,6,5,6,6,7,4,5,5,6,5,6,6,7,5,6,6,7,6,7,7,8};

static int Own_Correct_Equ_E_32f(Ipp32f num1, Ipp32f num2)
{
    if(((num1-num2)<=0.000001*IPP_MAX(num1,num2))&&((num1-num2)>=-0.000001*IPP_MAX(num1,num2)))
        return 1;
    else
        return 0;
}

Select_Gauss::Select_Gauss (int cdbk_size, bool is_row, int n_vec, float val, int num_clust, int spc_dim) {
   cdbksize      = cdbk_size;
   num_clusters  = num_clust;
   treshold      = val;
   num_vec       = n_vec;
   this->is_row  = is_row;
   space_dim     = spc_dim;
   space_dim4    = 0;
   gauss_num     = 0;
   max_len       = 0;
   states_num    = 0;
   count_all     = 0;
   count_one     = 0;
   pCdbk         = NULL;
   pIndx         = NULL;
   MixTable      = NULL;
   local_weights = NULL;
   pShortlist    = NULL;
}

Select_Gauss::~Select_Gauss(void) {
//   int j;

   if (pCdbk)         ippsCdbkFree_32f(pCdbk);
   if (pIndx)         ippsFree(pIndx);
   if (MixTable)      ippsFree(MixTable);
   if (local_weights) ippsFree(local_weights);
   if (pShortlist)    ippsFree(pShortlist);
   pCdbk         = NULL;
   pIndx         = NULL;
   MixTable      = NULL;
   local_weights = NULL;
   pShortlist    = NULL;
   cdbksize      = 0;
   num_clusters  = 0;
   treshold      = 0;
   num_vec       = 0;
   is_row        = 0;
   space_dim     = 0;
   space_dim4    = 0;
   gauss_num     = 0;
   max_len       = 0;
   states_num    = 0;
   count_all     = 0;
   count_one     = 0;
}

int Select_Gauss::Init_SelectGauss_Simple(int m_len, int m_gauss, int gss_num, float *tmpcodebook,
                                          float *weights, float *meanSet, Mixture *st, int st_num) {
   int i;

   if (cdbksize<2)              return -1;
   if (space_dim<=0)            return -1;

   max_len=m_len,
   max_gauss=m_gauss,
   gauss_num=gss_num;
   states=st;
   states_num=st_num;
   space_dim4=(space_dim+3)&(~3);

   if (num_clusters<1)          num_clusters=1;
   if (num_clusters>cdbksize)   num_clusters=cdbksize;
   if (num_vec>max_gauss)       num_vec=max_gauss;
   if (num_vec<0)               num_vec=0;
   if (treshold<1.0f)           treshold=1.0f;

   if (!(local_weights=ippsMalloc_32f(space_dim4))) {
      return -2;
   }
   ippsCopy_32f(weights,local_weights,space_dim4);


   if (ippsCdbkInitAlloc_WgtL2_32f(&pCdbk,tmpcodebook,weights,space_dim4,space_dim4,cdbksize,cdbksize,
                                    IPP_CDBK_FULL)!=ippStsNoErr) {
      return -4;
   }

   i=CommonInit(meanSet);

   if (i!=0) return -4;

   return 0;
}

int Select_Gauss::Init_SelectGauss_Full(int m_len, int m_gauss, int gss_num, float *tmpcodebook,
                                        float *t_weights, float *meanSet, Mixture *st,int st_num) {
   int i,j;
   float *g_means,*g_vars;

   if (cdbksize<2)              return -1;
   if (space_dim<=0)            return -1;

   max_len=m_len,
   max_gauss=m_gauss,
   gauss_num=gss_num;
   states=st;
   states_num=st_num;
   space_dim4=(space_dim+3)&(~3);

   if (num_clusters<1)          num_clusters=1;
   if (num_clusters>cdbksize)   num_clusters=cdbksize;
   if (num_vec>max_gauss)       num_vec=max_gauss;
   if (num_vec<0)               num_vec=0;
   if (treshold<1.0f)           treshold=1.0f;

   if (!(local_weights=ippsMalloc_32f(space_dim4*3))) {
      return -2;
   }

   g_means=local_weights+space_dim4;
   g_vars=g_means+space_dim4;

   if (ippsMeanVarColumn_32f_D2(meanSet,gauss_num,space_dim4,g_means,g_vars,space_dim4)!=ippStsNoErr) {
      return -3;
   }
   ippsSet_32f(1.0f,g_vars+space_dim,space_dim4-space_dim);

   ippsSet_32f(1.0f,t_weights,space_dim4);
   i=ippsDiv_32f_I(g_vars,t_weights,space_dim4);

   if(i!=0){                           // It means problem occured in ippsDiv_32f_I function
       if(i==6){                       // Zero value(s) of the divisor in the function Div
           i=0;
           for(j=0;j<space_dim4;j++){  // Should manually divide vector
               if(Own_Correct_Equ_E_32f(g_vars[j],0.0f)){ // g_vars[j]==0.0f
                   t_weights[j]=0.0f;
                   i++;
               }else{
                   t_weights[j]=1.0f/g_vars[j];
               }
           }

           if(i==space_dim4){          // All g_vars is equal ZERO so creating of codebook is impossible
               return -5;              // Only one cluster can be created so this example will be senseless
           }
       }else{
           return -6;
       }
   }


   if (ippsCdbkInitAlloc_WgtL2_32f(&pCdbk,meanSet,t_weights,space_dim4,space_dim4,gauss_num,cdbksize,
      IPP_CDBK_KMEANS_NUM)!=ippStsNoErr) {
      return -4;
   }

   ippsCopy_32f(g_means,local_weights,space_dim4);
   if (ippsGetCodebook_32f(pCdbk,tmpcodebook,space_dim4)!=0) {
      return -3;
   }

   i=CommonInit(meanSet);

   if (i!=0) return -4;

   return 0;
}

int Select_Gauss::CreateShortList(float *meanSet) {
   Ipp32s* tmpIndx;
   int i,j,k,l;
   CdbkState_32f* tmpCdbk;
   Ipp32f* tmpVec;
   Ipp32s* pNums;

   Listlen=((gauss_num+7)&(~7))/8;
   if (!(pShortlist=(Ipp8u**)ippsMalloc_8u((sizeof(Ipp8u*)+sizeof(Ipp8u)*Listlen)*cdbksize))) {
       return -2;
   }

   for (i=0; i<cdbksize; i++)
       pShortlist[i]=((Ipp8u*)(pShortlist+cdbksize))+i*Listlen;
   ippsSet_8u(0,pShortlist[0],cdbksize*Listlen);

   if (is_row) { //******* Row fill *******

       if (!(tmpVec=ippsMalloc_32f((space_dim4+max_gauss+1)*cdbksize)))
            return -2;
       tmpIndx=(Ipp32s*)(tmpVec+space_dim4*cdbksize);
       pNums=tmpIndx+max_gauss*cdbksize;

       if (ippsGetCodebook_32f(pCdbk,tmpVec,space_dim4)!=ippStsNoErr)
            return -4;

       if (num_vec>0) {

          if(num_vec>max_gauss) num_vec=max_gauss;
          ippsSet_32s(0,tmpIndx,cdbksize*num_vec);

          for (j=0; j<states_num; j++) {

              if (ippsCdbkInitAlloc_WgtL2_32f(&tmpCdbk,meanSet+space_dim4*states[j].compInd,
                 local_weights,space_dim4,space_dim4,states[j].compNum,states[j].compNum,
                 IPP_CDBK_FULL)!=ippStsNoErr) {
                 return -4;
              }
              if (num_vec>states[j].compNum) l=states[j].compNum;
              else                           l=num_vec;

              for (i=0; i<cdbksize; i++) {
                  if (ippsVQSingle_Sort_32f(tmpVec+space_dim4*i,tmpIndx+l*i,tmpCdbk,l)!=ippStsNoErr) {
                     return -4;
                  }
              }
              if (ippsFillShortlist_Row_1u(tmpIndx,states[j].compNum,l,pShortlist,cdbksize,Listlen,
                 states[j].compInd)!=ippStsNoErr) {
                 return -4;
              }
              ippsCdbkFree_32f(tmpCdbk);
          }

       } else {
          ippsSet_32s(0,tmpIndx,cdbksize*max_gauss);

          for (j=0; j<states_num; j++) {
              if (ippsCdbkInitAlloc_WgtL2_32f(&tmpCdbk,meanSet+space_dim4*states[j].compInd,
                 local_weights,space_dim4,space_dim4,states[j].compNum,states[j].compNum,
                 IPP_CDBK_FULL)!=ippStsNoErr) {
                 return -4;
              }
              for (i=0,k=0; i<cdbksize; i++) {
                 ippsVQSingle_Thresh_32f(tmpVec+space_dim4*i,tmpIndx+k,tmpCdbk,treshold,&pNums[i]);
                 k+=pNums[i];
              }
              ippsFillShortlist_RowVar_1u(tmpIndx,pNums,states[j].compNum,pShortlist,cdbksize,Listlen,
                 states[j].compInd);
              ippsCdbkFree_32f(tmpCdbk);
            }
       }
       if (tmpVec) ippsFree(tmpVec);

   } else {  //******* Column fill *******

       if (!(tmpIndx=ippsMalloc_32s(cdbksize)))
            return -2;

       if(num_vec>cdbksize) num_vec=cdbksize;

       for (i=0; i<gauss_num; i++) {
           if (num_vec>0) {
               ippsVQSingle_Sort_32f(meanSet+space_dim4*i,tmpIndx,pCdbk,num_vec);
               j=num_vec;
           } else {
               ippsVQSingle_Thresh_32f(meanSet+space_dim4*i,tmpIndx,pCdbk,treshold,&j);
           }
           ippsFillShortlist_Column_1u(tmpIndx,j,pShortlist,cdbksize,Listlen,i,1);
       }
       if (tmpIndx) ippsFree(tmpIndx);
   }

   return 0;
}

int Select_Gauss::CommonInit(float *means) {

   if (!(pIndx=(int*)ippsMalloc_32s(max_len*num_clusters))) {
      return -2;
   }

   if (!(MixTable=(Ipp8u*)ippsMalloc_8u(max_len*max_gauss))) {
      return -2;
   }

   ippsSet_8u(0x00,MixTable,max_len*max_gauss);

    return CreateShortList(means);
}

int Select_Gauss::WriteVQ(float* feature, int shift, int len) {
   int i;

   for (i=0; i<len; i++) {
      if (ippsVQSingle_Sort_32f(feature+space_dim4*i, pIndx+(shift+i)*num_clusters, pCdbk, num_clusters)!=ippStsNoErr)
         return -1;
    }

    return 0;
}


int Select_Gauss::ShiftIndxBegin(int len) {
   int i,j;

   for (i=0; i<len; i++)
      for (j=0; j<num_clusters; j++)
         pIndx[i*num_clusters+j]=pIndx[max_len*num_clusters-len*num_clusters+i*num_clusters+j];
    return 0;
}

int Select_Gauss::FillSign(int esnmix, int shiftbits, int num,int shift) {
    int i,k;

    ippsBuildSignTable_8u1u(pIndx+shift*num_clusters, num_clusters, (const Ipp8u**)pShortlist, cdbksize,
      Listlen, shiftbits, MixTable, num, esnmix);

    k=(esnmix+7)/8;
    count_all+=num*esnmix;
    for(i=0;i<num*k;i++)
        count_one+=CountArray[MixTable[i]];
    return 0;
}