home *** CD-ROM | disk | FTP | other *** search
/ ftp.disi.unige.it / 2015-02-11.ftp.disi.unige.it.tar / ftp.disi.unige.it / pub / .person / BarlaA / sw / OLD / Simo / SVM_mono / kernel_old.c < prev    next >
C/C++ Source or Header  |  2002-06-25  |  6KB  |  267 lines
  1. #include <math.h>
  2. #include <stdio.h>
  3. #include <string.h>
  4. #include "nrutil.h"
  5. #include "patt_rec.h"
  6.  
  7. extern int    dim;
  8. extern int    degree;
  9. extern double    sigma, normalization;
  10. extern int      dim_r,dim_c, delta_g, delta_r, delta_c;
  11.  
  12. double linear_kernel(double *v1, double *v2)
  13. {
  14.     register  i;
  15.     double aux = 0.;
  16.  
  17.     for (i = dim; i > 0; i--)
  18.         aux += *v1++ * *v2++;
  19. /*    fprintf(stderr, "%lf\n", aux/normalization);*/
  20.     return aux/normalization;
  21. }
  22.  
  23. double negative_kernel(double *v1, double *v2)
  24. {
  25.     register  i;
  26.     double aux = 0., tmp;
  27.  
  28.     
  29.     for (i = dim; i > 0; i--) {
  30.         tmp = *v1++ - *v2++;
  31.             aux -= tmp * tmp;
  32.       }
  33.     
  34.     /*
  35.     for (i=dim; i>0; i--) {
  36.       tmp = - (*v1)*(*v1) - (*v2)*(*v2) + 2*(*v1++)(*v2++);
  37.       aux += tmp;
  38.     }
  39.     */
  40.  
  41.  
  42.     return (aux/normalization);
  43. }
  44.  
  45. double poly_kernel(double *v1, double *v2)
  46. {
  47.     register  i;
  48.     double aux = 0.;
  49.     
  50.     for (i = dim; i > 0; i--)
  51.         aux += *v1++ * *v2++;
  52. /*    fprintf(stderr, "%lf\n", pow((aux/normalization + 1.),(double)degree)-1.);*/
  53.     return(pow((aux/normalization + 1.),(double)degree)-1.);
  54. }
  55.  
  56. double rbf_kernel(double *v1, double *v2)
  57. {
  58.     register  i;
  59.     double aux, res = 0.;
  60.     
  61.     for (i = dim; i > 0; i--)
  62.     {
  63.         aux = *v1++ - *v2++;
  64.         res += aux * aux;
  65.     }
  66. /*    fprintf(stderr, "%lf\n", exp(-res/(2*sigma*sigma)));*/
  67.     return exp(-res / (2*sigma*sigma));
  68. }
  69.  
  70. double mv_kernel(double *mod,double *imm)
  71. {
  72.   register i;
  73.   double normal_imm = 0.; 
  74.   double normal_mod = 0.; 
  75.   double mod2[dim],imm2[dim];    /*vettori ausiliarii*/
  76.   double inter = 0.;
  77.   
  78.   for(i=0; i<dim; i++){ 
  79.     //fprintf(stdout,"uffa %lf %lf \n",imm[i],mod[i]);
  80.     normal_imm += imm[i];
  81.     normal_mod += mod[i];
  82.   }
  83.   
  84.   //fprintf(stdout,"%lf %lf\n",normal_imm,normal_mod);
  85.   for(i=0; i<dim; i++){
  86.     imm2[i] = imm[i]/normal_imm;
  87.     mod2[i] = mod[i]/normal_mod;
  88.     inter +=  ilmiomin(imm2[i],mod2[i]);
  89.   }
  90.   
  91.   
  92.   return(inter); 
  93.  
  94.   
  95. } 
  96.  
  97. double haus_kernel(double *v1, double *v2)
  98. {
  99.    register int big_g;
  100.   int i,j,k,ii,jj,kk;
  101.   int max_val, counter, tmp_counter;
  102.   int st_k;
  103.   int check=0;
  104.  
  105.   int verbose, occlusion, max_noise;
  106.  
  107.   int type;
  108.   int w_r, w_c;
  109.   double *data_i, *model_data;
  110.   double *data_s, *base, *startp;
  111.   byte pop[ 256 ], *grey;
  112.  
  113.   double  *image_s;
  114.   long *addr,*next;
  115.   int **surface, **row_surface, **col_surface, **and,*buffer;
  116.  
  117.   long tot=0;
  118.   int start_h;
  119.   int stop_h;
  120.   int start_k;
  121.   int stop_k;
  122.  
  123.   model_data= v1; //model first
  124.   data_i = v2;    //image second
  125.  
  126.   w_r= dim_r;
  127.   w_c = dim_c;
  128.  
  129.   /* BEGIN initialization                 */
  130.  
  131.   for (i = 0; i < 256; i++)
  132.     pop[i] = 0;
  133.  
  134.   /*      pop[ ind ] is 1 if the grey value ind is present in the model     */
  135.  
  136.   for ( j = 0; j < w_r; j++ )
  137.     for ( k = 0; k < w_c; k++ ) 
  138.       {
  139.       pop[((byte) model_data[ j * dim_c + k ]) ] = 1;
  140.  
  141.     }
  142.  
  143.   /*  big_g will contain the number  of grey levels of the model          */
  144.  
  145.   big_g = 0;
  146.  
  147.   for (i = 0; i < 256; i++)
  148.     if (pop[ i ])
  149.       big_g++;
  150.  
  151.   grey = cvector(0, big_g - 1);
  152.  
  153.   j = 0;
  154.  
  155.   for (i = 0; i < 256; i++ )
  156.     if (pop[ i ] ) {
  157.       grey[ j ] = i;
  158.       pop[ i ]=j;
  159.       j++;
  160.     }
  161.  
  162.   /* addr is a buffer with as many elements as the model. Used to optimize the access to
  163.      the right element of the image data structure                                     */
  164.  
  165.   addr = lvector( 0, w_r*w_c );
  166.   for( i=0; i < w_r; i++ )
  167.     for( j=0; j < w_c; j++ ) {
  168.     addr[ i*w_c+j ] = pop[((byte)model_data[i*dim_c+j])]*dim_c*dim_r+i*dim_c+j;
  169.     }
  170.   addr[ w_c*w_r ]= -2;
  171.  
  172.  
  173.   data_s = dvector(0, big_g*dim_r*dim_c-1);
  174.   for(i=0;i<big_g*dim_r*dim_c-1;i++)
  175.     data_s[i]=0;
  176.  
  177. /* occhio allo sfondo = 180 */
  178.   for (j = 0; j < dim_r ; j++)
  179.     for (k = 0; k < dim_c; k++)
  180.        for (i = ((byte)data_i[j*dim_c+k]) -delta_g ; i <= ((byte)data_i[j*dim_c+k]) +delta_g ; i++)
  181.       for (jj = j - delta_r; jj <= j + delta_r; jj++)
  182.         for (kk = k - delta_c; kk <= k + delta_c; kk++) {
  183.           if ((i>=0)&&(i<256) &&(jj>=0)&&(jj<dim_r)&&(kk>=0)&&(kk<dim_c))
  184.         if (pop[i]!=0)
  185.           *(data_s+pop[i]*dim_c*dim_r+jj*dim_c+kk) = 1;
  186.       }
  187.  
  188.  
  189.   startp=data_s;
  190.   max_val = 0;
  191.  
  192.   counter=0;
  193.   next=addr;
  194.   base = startp;  //this command is equivalent to translating the model
  195.                               //along the image
  196.   for(;*next>=-1;) {
  197.  
  198.     if (*next!=-1) {
  199.       if (*(base+(*next++)))    //Jump in the right position of the image 3D structure
  200.                             //accessing the offsets stored in addr
  201.     counter++;              //Count the ones
  202.     }
  203.     else *next++;
  204.   }
  205.    
  206.  
  207.    free_dvector(data_s,0, big_g*dim_r*dim_c-1);
  208.    free_cvector(grey,0, big_g - 1);
  209.    free_lvector(addr, 0, w_r*w_c );
  210.  
  211.    /* AGGIUNTA SCOPIAZZAMENTO TIPI PROGETTO */
  212.    /* x e y sono le ccordinate del punto nell'immagine (i,j) */
  213.    
  214.    for(i=0;i<w_r;i++)//modello
  215.      for(j=0;j<w_c;j++){//immagine
  216.        start_h = (x-delta_c)+i;
  217.        stop_h  = (x+delta_c)+i;
  218.        start_k = (y-delta_r)+j;
  219.        stop_k  = (y+delta_r)+j;
  220.        for(h=start_h;h<stop_h;h++)
  221.      for(k=start_k;k<stop_k;k++){
  222.        
  223.        if ( abs((ActualImage->m_Matrix[h][k])-(ActualModel->m_Matrix[i][j])) <= epsilon ){
  224.          tot++;
  225.          h=stop_h;
  226.          k=stop_k;
  227.        }
  228.      }
  229.      }
  230.    return tot/(float)(m_w_r*m_w_c);
  231.    return((double)counter/normalization);
  232. }
  233.  
  234.  
  235.  
  236. static double (*kernel_v[])() = {linear_kernel, negative_kernel,poly_kernel, rbf_kernel, mv_kernel, haus_kernel};
  237.  
  238. double kernel(int k, double *v1, double *v2)
  239. {
  240.     switch (k)
  241.     {
  242.         case LINEAR_KERNEL:
  243.             return linear_kernel(v1, v2);
  244.         case NEGATIVE_KERNEL:
  245.             return negative_kernel(v1, v2);
  246.         case POLY_KERNEL:
  247.             return poly_kernel(v1, v2);
  248.         case RBF_KERNEL:
  249.             return rbf_kernel(v1, v2);
  250.         case MV_KERNEL:
  251.             return mv_kernel(v1, v2);
  252.         case HAUS_KERNEL:
  253.           {
  254.             return haus_kernel(v1, v2);
  255.           }
  256.                 case HAUS_KERNEL_SIMM:
  257.                   {
  258.                     return (haus_kernel(v1, v2)+haus_kernel(v2, v1))/2;
  259.                    }
  260.                 fprintf(stderr, "Unknown kernel\n");
  261.             exit(-1);
  262.     }
  263.     return (*kernel_v[k])(v1, v2);
  264. }
  265.  
  266.  
  267.