levelines_kernels.cu

   1
   2 // chemins des lignes de niveaux
   3 // longueur = 4 pixels
   4 // une ligne = un chemin
   5
   6 __constant__ int pathDi[PSIZE_I][PSIZE_J-1] =
   7   {
   8     // Q1
   9         {  -1, -1, -1},       // 90
  10         {  -1, -1, -1},       // 75
  11         {  -1, -1, -1},       // 60
  12         {  -1, -1, -1},       // 45
  13         {  -1,  0, -1},       // 30
  14         {   0, -1,  0},       // 15
  15         // Q4
  16         {   0,  0,  0},       // 0
  17         {   0,  1,  0},       // 345
  18         {   1,  0,  1},       // 330
  19         {   1,  1,  1},       // 315
  20         {   1,  1,  1},       // 300
  21         {   1,  1,  1},       // 285
  22         // Q3
  23         {   1,  1,  1},       // 270
  24         {   1,  1,  1},       // 255
  25         {   1,  1,  1},       // 240
  26         {   1,  1,  1},       // 225
  27         {   1,  0,  1},       // 210
  28         {   0,  1,  0},       // 195
  29         // Q2
  30         {   0,  0,  0},       // 180
  31         {   0, -1,  0},       // 165
  32         {  -1,  0, -1},       // 150
  33         {  -1, -1, -1},       // 135
  34         {  -1, -1, -1},       // 120
  35         {  -1, -1, -1}        // 105
  36   } ;     //
  37
  38 __constant__ int pathDj[PSIZE_I][PSIZE_J-1] =
  39   {
  40         // Q1
  41         {  0,  0,  0},       // 90
  42         {  0,  1,  0},       // 75
  43         {  1,  0,  1},       // 60
  44         {  1,  1,  1},       // 45
  45         {  1,  1,  1},       // 30
  46         {  1,  1,  1},       // 15
  47         // Q4
  48         {  1,  1,  1},       // 0
  49         {  1,  1,  1},       // 345
  50         {  1,  1,  1},       // 330
  51         {  1,  1,  1},       // 315
  52         {  1,  0,  1},       // 300
  53         {  0,  1,  0},       // 285
  54         // Q3
  55         {  0,  0,  0},       // 270
  56         {  0, -1,  0},       // 255
  57         { -1,  0, -1},       // 240
  58         { -1, -1, -1},       // 225
  59         { -1, -1, -1},       // 210
  60         { -1, -1, -1},       // 195
  61         // Q2
  62         { -1, -1, -1},       // 180
  63         { -1, -1, -1},       // 165
  64         { -1, -1, -1},       // 150
  65         { -1, -1, -1},       // 135
  66         { -1,  0, -1},       // 120
  67         {  0, -1,  0}        // 105
  68   } ;
  69
  70
  71 // declare texture reference for 2D int texture
  72 texture<int, 2, cudaReadModeElementType> tex_img_in ;
  73 texture<int, 2, cudaReadModeElementType> tex_img_estim ;
  74 texture<int, 2, cudaReadModeElementType> tex_img_lniv ;
  75
  76
  77 __constant__ float tangente[] = {0.000, 0.268, 0.577, 1.000} ;
  78
  79 __global__ void kernel_calcul_paths( int2 * d_paths, unsigned int r){
  80
  81   unsigned int idpath = 0 ;
  82   int ic, jc, iprec, jprec ;
  83   float offset = 0.5 ;
  84   unsigned int basepath = 0 ;
  85
  86   // Q1 inf
  87   for (int a=0 ; a< 4 ; a++){        // les 4 angles 0,15,30 et 45
  88         for (int p=0 ; p< r ; p++){      // les r points
  89           ic = r-1 - floor(tangente[a]*p + offset) ;
  90           if ( p > 0 ){
  91                 d_paths[idpath*(r-1)+p-1].x = ic - iprec ;
  92                 d_paths[idpath*(r-1)+p-1].y = 1 ;
  93           }
  94           iprec = ic ;
  95         }
  96         idpath++ ;
  97   }
  98   // Q1 sup
  99   for (int a=2 ; a>0 ; a--){         // les 2 angles 60 et 75
 100         for (int p=0 ; p< r ; p++){      // les r points
 101           jc = floor(tangente[a]*p + offset) ;
 102           if ( p > 0 ){
 103                 d_paths[idpath*(r-1)+p-1].x = -1 ;
 104                 d_paths[idpath*(r-1)+p-1].y = jc - jprec ;
 105           }
 106           jprec = jc ;
 107         }
 108         idpath++ ;
 109   }
 110
 111   // Q2
 112   basepath += 6 ;
 113   for (int a=0 ; a< 6 ; a++){         // les 6 angles 90,105,120,135,150,165
 114         for (int p=0 ; p<r-1 ; p++){      // les r points
 115           d_paths[idpath*(r-1)+p].x = -d_paths[(idpath - basepath)*(r-1)+p].y ;
 116           d_paths[idpath*(r-1)+p].y =  d_paths[(idpath - basepath)*(r-1)+p].x ;
 117           }
 118         idpath++ ;
 119   }
 120
 121   // Q3
 122   basepath += 6 ;
 123   for (int a=0 ; a< 6 ; a++){         // les 6 angles 180,195,210,225,240,255
 124         for (int p=0 ; p<r-1 ; p++){      // les r points
 125           d_paths[idpath*(r-1)+p].x = -d_paths[(idpath - basepath)*(r-1)+p].x ;
 126           d_paths[idpath*(r-1)+p].y = -d_paths[(idpath - basepath)*(r-1)+p].y ;
 127           }
 128         idpath++ ;
 129   }
 130
 131   // Q4
 132   basepath += 6 ;
 133   for (int a=0 ; a< 6 ; a++){         // les 6 angles 270,285,300,315,330,345
 134         for (int p=0 ; p<r-1 ; p++){      // les r points
 135           d_paths[idpath*(r-1)+p].x =  d_paths[(idpath - basepath)*(r-1)+p].y ;
 136           d_paths[idpath*(r-1)+p].y = -d_paths[(idpath - basepath)*(r-1)+p].x ;
 137           }
 138         idpath++ ;
 139   }
 140 }
 141
 142
 143 __global__ void kernel_init_estim_from_img_in(unsigned int * d_estim, unsigned int L, unsigned int H, unsigned int r){
 144   // coordonnes du point dans l'image
 145   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 146   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 147   unsigned int pos = i*L +j ;
 148   unsigned int ic , jc, ng ;
 149
 150   if( (i>r)&&(i<H-r)&&(j>r)&&(j<L-r) ){
 151         ng = 0 ;
 152         for (ic = i - r ; ic <= i + r ; ic++){
 153           for (jc = j - r ; jc <= j + r ; jc++){
 154                 ng += tex2D(tex_img_in, jc, ic ) ;
 155           }
 156         }
 157         d_estim[ pos ] = ng/((2*r+1)*(2*r+1)) ;
 158   }
 159   /*
 160   else
 161   d_estim[ pos ] = tex2D(tex_img_in, j, i) ;
 162   */
 163 }
 164
 165 __global__ void kernel_neutre_estim_from_img_in( unsigned int * d_estim, unsigned int L, unsigned int H ){
 166   // coordonnes du point dans l'image
 167   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 168   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 169   unsigned int pos = i*L +j ;
 170
 171   d_estim[ pos ] = tex2D(tex_img_in, j, i) ;
 172
 173
 174 }
 175
 176 __global__ void kernel_init_estim_from_img_in_global_mem(unsigned int * d_data, unsigned int * d_estim,
 177                                                                                                                  unsigned int L, unsigned int H, unsigned int r){
 178   // coordonnes du point dans l'image
 179   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 180   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 181   unsigned int pos = i*L +j ;
 182   unsigned int ic , jc, ng ;
 183
 184   if( (i>r)&&(i<H-r)&&(j>r)&&(j<L-r) ){
 185         ng = 0 ;
 186         for (ic = i - r ; ic <= i + r ; ic++){
 187           for (jc = j - r ; jc <= j + r ; jc++){
 188                 ng += d_data[ ic*L + jc ] ;
 189           }
 190         }
 191         d_estim[ pos ] = ng/((2*r+1)*(2*r+1)) ;
 192   }
 193 }
 194
 195 __global__ void kernel_estim_next_step_global_mem(unsigned int * d_estim, unsigned int * d_lniv, unsigned int L, unsigned int H, unsigned int p){
 196   // coordonnes du point dans l'image
 197   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 198   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 199   unsigned int pos = i*L +j ;
 200
 201   d_estim[ pos ] = ( tex2D(tex_img_in, j, i) + p*d_lniv[ pos ] )/(1+p) ;
 202
 203 }
 204
 205 __global__ void kernel_estim_next_step_texture(unsigned int * d_estim, unsigned int L, unsigned int H, unsigned int p){
 206   // coordonnes du point dans l'image
 207   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 208   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 209   unsigned int pos = i*L +j ;
 210
 211   d_estim[ pos ] = ( tex2D(tex_img_in, j, i) + p*tex2D(tex_img_lniv, j, i ))/(1+p) ;
 212
 213 }
 214
 215 __global__ void kernel_estim_next_step_global_mem(unsigned int * d_estim, unsigned int * d_lniv, unsigned int * d_data,
 216                                                                                                   unsigned int L, unsigned int H, unsigned int p){
 217   // coordonnes du point dans l'image
 218   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 219   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 220   unsigned int pos = i*L +j ;
 221
 222   d_estim[ pos ] = ( d_data[ pos ] + p*d_lniv[ pos ])/(1+p) ;
 223
 224 }
 225
 226 __global__ void kernel_levelines_global_mem(unsigned int * img_in, unsigned int * img_out, unsigned int L, unsigned int H)
 227 {
 228   // coordonnes du point dans l'image
 229   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 230   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 231   //unsigned int spos = threadIdx.x * blockDim.y + threadIdx.y ;
 232
 233   // nb de points par chemin
 234   int lpath =  PSIZE_J ;
 235   unsigned int ic, jc, zc, pos = i*L+j;
 236   int idpath, idpix ;
 237   unsigned int mse_min, mse_cur, val ;
 238   uint2 mse ;
 239
 240
 241   if((i>lpath)&&(i<H-lpath)&&(j>lpath)&&(j<L-lpath)){
 242         for( idpath=0; idpath < PSIZE_I ; idpath++) {
 243           ic = i ;
 244           jc = j ;
 245           pos = ic*L + jc ;
 246           zc = img_in[ pos ] ;
 247           mse.x = zc ;
 248           mse.y = zc*zc ;
 249           for( idpix=0; idpix < lpath-1 ; idpix++ ) {
 250                 ic += pathDi[idpath][idpix] ;
 251                 jc += pathDj[idpath][idpix] ;
 252                 pos = ic*L + jc ;
 253                 zc = img_in[ pos ] ;
 254                 mse.x += zc ;
 255                 mse.y += zc*zc ;
 256           }
 257           // critere de selection du chemin ( SUM_(X2) - SUM_(X)2 / lpath )
 258           // a ameliorer pour vitesse
 259           mse_cur = ( mse.y - ( mse.x / lpath ) * mse.x ) ;
 260           if (idpath == 0) {
 261                 mse_min = mse_cur ;
 262                 val = mse.x ;
 263           } else {
 264                 if ( mse_cur < mse_min )  {
 265                   mse_min = mse_cur ;
 266                   val = mse.x ;
 267                 }
 268           }
 269         }
 270         img_out[ i*L + j ] = val / lpath ;
 271   }
 272 }
 273
 274 __global__ void kernel_levelines_texture(unsigned int * img_out, unsigned int L, unsigned int H)
 275 {
 276   // coordonnes du point dans l'image
 277   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 278   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 279   //unsigned int spos = threadIdx.x * blockDim.y + threadIdx.y ;
 280
 281   // nb de points par chemin
 282   int lpath =  PSIZE_J ;
 283   unsigned int ic, jc, zc ;
 284   int idpath, idpix ;
 285   unsigned int mse_min, mse_cur, val ;
 286   uint2 mse ;
 287
 288
 289   if( (i >=lpath )&&( i<= (H-lpath) )&&( j>=lpath )&&( j<=(L-lpath) ) ){
 290         for( idpath=0; idpath < PSIZE_I ; idpath++) {
 291           ic = i ;
 292           jc = j ;
 293           zc = tex2D(tex_img_estim, j, i) ;
 294           mse.x = zc ;
 295           mse.y = zc*zc ;
 296           for( idpix=0; idpix < lpath-1 ; idpix++ ) {
 297                 ic += pathDi[idpath][idpix] ;
 298                 jc += pathDj[idpath][idpix] ;
 299                 zc = tex2D(tex_img_estim, jc, ic) ;
 300                 mse.x += zc ;
 301                 mse.y += zc*zc ;
 302           }
 303           // critere de selection du chemin ( SUM_(X2) - SUM_(X)2 / lpath )
 304           // a ameliorer pour vitesse
 305           mse_cur = ( mse.y - ( mse.x*mse.x / lpath ) ) ;
 306           if (idpath == 0) {
 307                 mse_min = mse_cur ;
 308                 val = mse.x ;
 309           } else {
 310                 if ( mse_cur < mse_min )  {
 311                   mse_min = mse_cur ;
 312                   val = mse.x ;
 313                 }
 314           }
 315         }
 316         img_out[ i*L + j ] = val / lpath ;
 317   }
 318   else img_out[ i*L + j ] = 0 ;
 319 }
 320
 321
 322 __global__ void kernel_dir_levelines_texture(unsigned int * dir, unsigned int L, unsigned int H)
 323 {
 324   // coordonnes du point dans l'image
 325   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 326   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;
 327   //unsigned int spos = threadIdx.x * blockDim.y + threadIdx.y ;
 328
 329   // nb de points par chemin
 330   int lpath =  PSIZE_J ;
 331   unsigned int ic, jc, zc ;
 332   int idpath, idpix ;
 333   unsigned int mse_min, mse_cur, val ;
 334   uint2 mse ;
 335
 336
 337   if( (i >=lpath )&&( i<= (H-lpath) )&&( j>=lpath )&&( j<=(L-lpath) ) ){
 338         for( idpath=0; idpath < PSIZE_I ; idpath++) {
 339           ic = i ;
 340           jc = j ;
 341           zc = tex2D(tex_img_estim, j, i) ;
 342           mse.x = zc ;
 343           mse.y = zc*zc ;
 344           for( idpix=0; idpix < lpath-1 ; idpix++ ) {
 345                 ic += pathDi[idpath][idpix] ;
 346                 jc += pathDj[idpath][idpix] ;
 347                 zc = tex2D(tex_img_estim, jc, ic) ;
 348                 mse.x += zc ;
 349                 mse.y += zc*zc ;
 350           }
 351           // critere de selection du chemin ( SUM_(X2) - SUM_(X)2 / lpath )
 352           // a ameliorer pour vitesse
 353           mse_cur = ( mse.y - ( mse.x*mse.x / lpath ) ) ;
 354           if (idpath == 0) {
 355                 mse_min = mse_cur ;
 356                 val = idpath ;
 357           } else {
 358                 if ( mse_cur < mse_min )  {
 359                   mse_min = mse_cur ;
 360                   val = idpath ;
 361                 }
 362           }
 363         }
 364         dir[ i*L + j ] = val ;
 365   }
 366
 367 }
 368
 369
 370
 371
 372 __global__ void kernel_trace_levelines_texture(unsigned int * dir, unsigned int * img_out,
 373                                                                            unsigned int L, unsigned int H, unsigned int pas, unsigned int ng){
 374   // coordonnes du point dans l'image
 375   unsigned int i = blockIdx.x*blockDim.x + threadIdx.x;
 376   unsigned int j = blockIdx.y*blockDim.y + threadIdx.y;;
 377
 378   // nb de points par chemin
 379   int lpath =  PSIZE_J ;
 380   unsigned int ic, jc, idpix ;
 381   unsigned int idpath  ;
 382
 383   img_out[ i*L + j ] = tex2D(tex_img_estim, j, i ) ;
 384
 385   if ( !(i%pas+j%pas)&&(i>=lpath)&&(i<=H-lpath)&&(j>=lpath)&&(j<=L-lpath) ){
 386         ic = i ;
 387         jc = j ;
 388         idpath = dir[ic*L+jc] ;
 389         img_out[ ic*L+jc ] = ng ;
 390         for ( idpix=0 ; idpix < lpath-1 ; idpix++ ){
 391           ic += pathDi[idpath][idpix] ;
 392           jc += pathDj[idpath][idpix] ;
 393           img_out[ ic*L + jc ] = ng ;
 394           }
 395   }
 396
 397 }