Go to the source code of this file.
Defines | |
| #define | _LAYER_ |
Functions | |
| void | process_layer (float **noisy, int Ni, int Nj, float **denoised, int pi, int pj, int di, int dj, float threshold) |
| void | process_layer_w_signif_sets (float **noisy, int Ni, int Nj, float **denoised, int pi, int pj, int di, int dj, float threshold) |
|
|
|
|
||||||||||||||||||||||||||||||||||||||||
|
Definition at line 15 of file layer.c. Referenced by fill_layers.
00017 {
00018 int i,j,k,l,cnt;
00019
00020 int **normalization;
00021 float **buffer_t,**buffer_i,**H,**buffer;
00022
00023 int dct_copy_sti,dct_copy_endi;
00024 int dct_copy_stj,dct_copy_endj;
00025 int offset=DCT_SIZE-1;
00026 int ovrlp_factor=(DCT_SIZE+1)/2;
00027 int starti,startj;
00028 int oi,oj;
00029
00030 buffer_t=allocate_2d_float(DCT_SIZE,DCT_SIZE,1);
00031 buffer_i=allocate_2d_float(DCT_SIZE,DCT_SIZE,1);
00032
00033 // We will be doing fully translation invariant (DCT_SIZE x DCT_SIZE) DCTs.
00034 // At image broders not all DCT_SIZE^2 DCTs contribute, so we need
00035 // a mechanism of counting the contribution. The way I am handling
00036 // this is a very general/redundant way of doing it. But I use it in
00037 // other ideas/code as well (see my "Weighted Overcomplete Denoising" paper).
00038 //
00039 // counts the number of transforms that contributed to denoising at every pixel.
00040 normalization=allocate_2d_int(di,dj,1);
00041
00042 // accumulates the denoised results at every pixel for final normalization.
00043 buffer=allocate_2d_float(di,dj,1);
00044
00045 H=gen_1d_dct(DCT_SIZE);
00046
00047 starti=pi-offset;
00048 startj=pj-offset;
00049
00050 // traverse all (DCT_SIZE x DCT_SIZE) blocks overlapping the missing region
00051 // for overcomplete denoising.
00052 for(i=starti;i<pi+di;i++) {
00053
00054 // check out of bounds.
00055 if((i<0)||(i+DCT_SIZE>Ni))
00056 continue;
00057
00058 //calculate overlap with missing region, see journal paper.
00059 if(i<pi) {
00060
00061 oi=i+DCT_SIZE-pi;
00062 dct_copy_sti=DCT_SIZE-oi;
00063 }
00064 else {
00065
00066 oi=pi+di-i;
00067 dct_copy_sti=0;
00068 }
00069
00070 if(oi>DCT_SIZE)
00071 oi=DCT_SIZE;
00072
00073 if(oi>di)
00074 oi=di;
00075
00076 dct_copy_endi=dct_copy_sti+oi;
00077
00078 for(j=startj;j<pj+dj;j++) {
00079
00080 // check out of bounds.
00081 if((j<0)||(j+DCT_SIZE>Nj))
00082 continue;
00083
00084 //calculate overlap with missing region.
00085 if(j<pj) {
00086
00087 oj=j+DCT_SIZE-pj;
00088 dct_copy_stj=DCT_SIZE-oj;
00089 }
00090 else {
00091
00092 oj=pj+dj-j;
00093 dct_copy_stj=0;
00094 }
00095
00096 if(oj>DCT_SIZE)
00097 oj=DCT_SIZE;
00098
00099 if(oj>dj)
00100 oj=dj;
00101
00102 dct_copy_endj=dct_copy_stj+oj;
00103
00104 // check amount of overlap, see journal paper.
00105 if((oi>ovrlp_factor)&&(oj>ovrlp_factor))
00106 continue;
00107
00108 // one block forward trf
00109 trf_one_block_dct_2d(noisy,i,j,buffer_t,0,0,DCT_SIZE,H);
00110
00111 // hard threshold
00112 cnt=hard_threshold(buffer_t,DCT_SIZE,DCT_SIZE,threshold);
00113
00114 // inverse trf
00115 inv_trf_one_block_dct_2d(buffer_t,0,0,buffer_i,0,0,DCT_SIZE,H);
00116
00117 // accumulate results for overcomplete denoising.
00118 for(k=dct_copy_sti;k<dct_copy_endi;k++) {
00119
00120 for(l=dct_copy_stj;l<dct_copy_endj;l++) {
00121
00122 normalization[i-pi+k][j-pj+l]++;
00123 buffer[i-pi+k][j-pj+l]+=buffer_i[k][l];
00124 }
00125 }
00126 }
00127 }
00128
00129 // left/right portion of the single pixel thick layer.
00130 for(i=pi;i<pi+di;i++) {
00131
00132 for(j=pj;j<pj+dj;j+=dj-1) {
00133
00134 // normalize accumulated results and copy into place.
00135 if(normalization[i-pi][j-pj])
00136 buffer[i-pi][j-pj]/=normalization[i-pi][j-pj];
00137 else
00138 buffer[i-pi][j-pj]=noisy[i][j];
00139
00140 denoised[i][j]=buffer[i-pi][j-pj];
00141 }
00142 }
00143
00144 // top/bottom portion of the single pixel thick layer.
00145 for(i=pi;i<pi+di;i+=di-1) {
00146
00147 for(j=pj+1;j<pj+dj-1;j++) {
00148
00149 // normalize accumulated results and copy into place.
00150 if(normalization[i-pi][j-pj])
00151 buffer[i-pi][j-pj]/=normalization[i-pi][j-pj];
00152 else
00153 buffer[i-pi][j-pj]=noisy[i][j];
00154
00155 denoised[i][j]=buffer[i-pi][j-pj];
00156 }
00157 }
00158
00159 free_2d_float(buffer_t,DCT_SIZE);
00160 free_2d_float(buffer_i,DCT_SIZE);
00161 free_2d_float(H,DCT_SIZE);
00162 free_2d_int(normalization,di);
00163 free_2d_float(buffer,di);
00164
00165 }
|
|
||||||||||||||||||||||||||||||||||||||||
|
Definition at line 234 of file layer.c. Referenced by fill_layers.
00236 {
00237 int i,j,k,l;
00238
00239 int **normalization;
00240 float **buffer_t,**buffer_i,**H,**buffer;
00241
00242 int dct_copy_sti,dct_copy_endi;
00243 int dct_copy_stj,dct_copy_endj;
00244 int offset=DCT_SIZE-1;
00245 int ovrlp_factor=(DCT_SIZE+1)/2;
00246 int starti,startj;
00247 int oi,oj,si,sj;
00248 int my_Ni,my_Nj;
00249
00250 static int ****non_zero=NULL,kilroy=0;
00251 static float my_thres=0,**my_noisy=NULL;
00252 static int my_pi=-1,my_pj=-1,my_di=-1,my_dj=-1;
00253
00254 if(!kilroy) {
00255
00256 kilroy=1;
00257 non_zero=(int ****)malloc(DCT_SIZE*sizeof(int ***));
00258 for(i=0;i<DCT_SIZE;i++) {
00259
00260 non_zero[i]=(int ***)malloc(DCT_SIZE*sizeof(int **));
00261
00262 for(j=0;j<DCT_SIZE;j++) {
00263
00264 non_zero[i][j]=allocate_2d_int(Ni,Nj,0);
00265 }
00266 }
00267
00268 my_thres=threshold;
00269 my_pi=pi;
00270 my_pj=pj;
00271 my_di=di;
00272 my_dj=dj;
00273
00274 my_Ni=di+2*offset;
00275 my_Nj=dj+2*offset;
00276
00277 my_noisy=allocate_2d_float(Ni,Nj,0);
00278
00279 starti=pi-offset;
00280 startj=pj-offset;
00281 for(i=starti;i<pi+di+offset;i++) {
00282
00283 // check out of bounds.
00284 if((i<0)||(i>=Ni))
00285 continue;
00286
00287 for(j=startj;j<pj+dj+offset;j++) {
00288
00289 // check out of bounds.
00290 if((j<0)||(j>=Nj))
00291 continue;
00292
00293 my_noisy[i-starti][j-startj]=noisy[i][j];
00294 }
00295 }
00296
00297 determine_dct_significant_sets(my_noisy,my_Ni,my_Nj,non_zero,threshold);
00298 }
00299 else if((my_thres!=threshold)||(pi!=my_pi)||(pj!=my_pj)||(di!=my_di)||(dj!=my_dj)){
00300
00301 // uncomment this to see when significant sets change if INNER_ITER_COUNT > 1 in recover.c
00302 /* if((my_thres!=threshold))
00303 printf("process_layer_w_signif_sets: threshold changed, recomputing significant sets.\n");
00304 else
00305 printf("process_layer_w_signif_sets: layer changed, recomputing significant sets.\n");*/
00306
00307 my_thres=threshold;
00308 my_pi=pi;
00309 my_pj=pj;
00310 my_di=di;
00311 my_dj=dj;
00312
00313 my_Ni=di+2*offset;
00314 my_Nj=dj+2*offset;
00315
00316 starti=pi-offset;
00317 startj=pj-offset;
00318 for(i=starti;i<pi+di+offset;i++) {
00319
00320 // check out of bounds.
00321 if((i<0)||(i>=Ni))
00322 continue;
00323
00324 for(j=startj;j<pj+dj+offset;j++) {
00325
00326 // check out of bounds.
00327 if((j<0)||(j>=Nj))
00328 continue;
00329
00330 my_noisy[i-starti][j-startj]=noisy[i][j];
00331 }
00332 }
00333
00334 determine_dct_significant_sets(my_noisy,my_Ni,my_Nj,non_zero,threshold);
00335 }
00336
00337 buffer_t=allocate_2d_float(DCT_SIZE,DCT_SIZE,1);
00338 buffer_i=allocate_2d_float(DCT_SIZE,DCT_SIZE,1);
00339
00340 normalization=allocate_2d_int(di,dj,1);
00341 buffer=allocate_2d_float(di,dj,1);
00342
00343 H=gen_1d_dct(DCT_SIZE);
00344
00345
00346 starti=pi-offset;
00347 startj=pj-offset;
00348 for(i=starti;i<pi+di;i++) {
00349
00350 // check out of bounds.
00351 if((i<0)||(i+DCT_SIZE>Ni))
00352 continue;
00353
00354 //calculate overlap with missing region.
00355 if(i<pi) {
00356
00357 oi=i+DCT_SIZE-pi;
00358 dct_copy_sti=DCT_SIZE-oi;
00359 }
00360 else {
00361
00362 oi=pi+di-i;
00363 dct_copy_sti=0;
00364 }
00365
00366 if(oi>DCT_SIZE)
00367 oi=DCT_SIZE;
00368
00369 if(oi>di)
00370 oi=di;
00371
00372 dct_copy_endi=dct_copy_sti+oi;
00373
00374 for(j=startj;j<pj+dj;j++) {
00375
00376 // check out of bounds.
00377 if((j<0)||(j+DCT_SIZE>Nj))
00378 continue;
00379
00380 //calculate overlap with missing region.
00381 if(j<pj) {
00382
00383 oj=j+DCT_SIZE-pj;
00384 dct_copy_stj=DCT_SIZE-oj;
00385 }
00386 else {
00387
00388 oj=pj+dj-j;
00389 dct_copy_stj=0;
00390 }
00391
00392 if(oj>DCT_SIZE)
00393 oj=DCT_SIZE;
00394
00395 if(oj>dj)
00396 oj=dj;
00397
00398 dct_copy_endj=dct_copy_stj+oj;
00399
00400 // check amount of overlap.
00401 if((oi>ovrlp_factor)&&(oj>ovrlp_factor))
00402 continue;
00403
00404 // forward trf
00405 trf_one_block_dct_2d(noisy,i,j,buffer_t,0,0,DCT_SIZE,H);
00406
00407 si=(i-starti)%DCT_SIZE;
00408 sj=(j-startj)%DCT_SIZE;
00409
00410 for(k=0;k<DCT_SIZE;k++) {
00411
00412 for(l=0;l<DCT_SIZE;l++) {
00413
00414 // significant set already calculated, no thresholding required.
00415 if(!non_zero[si][sj][i+k-starti][j+l-startj]) {
00416
00417 buffer_t[k][l]=0;
00418 }
00419 }
00420 }
00421
00422 // inverse trf
00423 inv_trf_one_block_dct_2d(buffer_t,0,0,buffer_i,0,0,DCT_SIZE,H);
00424
00425 for(k=dct_copy_sti;k<dct_copy_endi;k++) {
00426
00427 for(l=dct_copy_stj;l<dct_copy_endj;l++) {
00428
00429 // accumulate
00430 normalization[i-pi+k][j-pj+l]++;
00431 buffer[i-pi+k][j-pj+l]+=buffer_i[k][l];
00432 }
00433 }
00434 }
00435 }
00436
00437 // left/right portion of the single pixel thick layer.
00438 for(i=pi;i<pi+di;i++) {
00439
00440 for(j=pj;j<pj+dj;j+=dj-1) {
00441
00442 if(normalization[i-pi][j-pj])
00443 buffer[i-pi][j-pj]/=normalization[i-pi][j-pj];
00444 else
00445 buffer[i-pi][j-pj]=noisy[i][j];
00446
00447 denoised[i][j]=buffer[i-pi][j-pj];
00448 }
00449 }
00450
00451 // top/bottom portion of the single pixel thick layer.
00452 for(i=pi;i<pi+di;i+=di-1) {
00453
00454 for(j=pj+1;j<pj+dj-1;j++) {
00455
00456 if(normalization[i-pi][j-pj])
00457 buffer[i-pi][j-pj]/=normalization[i-pi][j-pj];
00458 else
00459 buffer[i-pi][j-pj]=noisy[i][j];
00460
00461 denoised[i][j]=buffer[i-pi][j-pj];
00462 }
00463 }
00464
00465 free_2d_float(buffer_t,DCT_SIZE);
00466 free_2d_float(buffer_i,DCT_SIZE);
00467 free_2d_float(H,DCT_SIZE);
00468 free_2d_int(normalization,di);
00469 free_2d_float(buffer,di);
00470
00471 }
|
1.2.14 written by Dimitri van Heesch,
© 1997-2002