diff --git a/rtengine/CA_correct_RT.cc b/rtengine/CA_correct_RT.cc index c84598c99..60133eb0a 100644 --- a/rtengine/CA_correct_RT.cc +++ b/rtengine/CA_correct_RT.cc @@ -5,7 +5,7 @@ // copyright (c) 2008-2010 Emil Martinec // // -// code dated: June 2, 2010 +// code dated: June 14, 2010 // // CA_correct_RT.cc is free software: you can redistribute it and/or modify // it under the terms of the GNU General Public License as published by @@ -100,10 +100,8 @@ int RawImageSource::LinEqSolve(int c, int dir, int nDim, float* pfMatr, float* p void RawImageSource::CA_correct_RT() { #define TS 256 // Tile size -//#define polyord 4 // max order of fit monomial -//#define numpar 16 // number of fit parameters = SQR(polyord) -#define border 8 -#define border2 16 +//#define border 8 +//#define border2 16 #define PIX_SORT(a,b) { if ((a)>(b)) {temp=(a);(a)=(b);(b)=temp;} } #define SQR(x) ((x)*(x)) @@ -113,18 +111,20 @@ void RawImageSource::CA_correct_RT() { float (*Gtmp); Gtmp = (float (*)) calloc ((height)*(width), sizeof *Gtmp); - int polyord=2, numpar=4, numblox[3]={0,0,0}; + static const int border=8; + static const int border2=16; + int polyord=4, numpar=16, numblox[3]={0,0,0}; - //static const int border=8; int rrmin, rrmax, ccmin, ccmax; int top, bottom, left, right, row, col; int rr, cc, rr1, cc1, c, indx, indx1, i, j, k, m, n, dir; + int areawt[2][3]; int GRBdir[2][3], offset[2][3]; int shifthfloor[3], shiftvfloor[3], shifthceil[3], shiftvceil[3]; - int vblsz, hblsz, vblock, hblock; + int vblsz, hblsz, vblock, hblock, vz1, hz1; //int verbose=1; int res; - int v1=TS, v2=2*TS, v3=3*TS, v4=4*TS;//, p1=-TS+1, p2=-2*TS+2, p3=-3*TS+3, m1=TS+1, m2=2*TS+2, m3=3*TS+3; + static const int v1=TS, v2=2*TS, v3=3*TS, v4=4*TS;//, p1=-TS+1, p2=-2*TS+2, p3=-3*TS+3, m1=TS+1, m2=2*TS+2, m3=3*TS+3; float eps=1e-10; //tolerance to avoid dividing by zero @@ -132,13 +132,13 @@ void RawImageSource::CA_correct_RT() { float coeff[2][3][3], CAshift[2][3]; float polymat[3][2][1296], shiftmat[3][2][36], fitparams[3][2][36]; float shifthfrac[3], shiftvfrac[3], temp, p[9]; - float gdiff, deltgrb, denom, Ginthfloor, Ginthceil, Gint, gradwt; + float gdiff, deltgrb, denom, Ginthfloor, Ginthceil, Gint, RBint, gradwt; float grbdiffinthfloor, grbdiffinthceil, grbdiffint, grbdiffold; float blockave[2][3]={{0,0,0},{0,0,0}}, blocksqave[2][3]={{0,0,0},{0,0,0}}, blockdenom[2][3]={{0,0,0},{0,0,0}}, blockvar[2][3]; - + float glpfh, glpfv, ghpfh, ghpfv; - static const float gaussg[5] = {0.171582, 0.15839, 0.124594, 0.083518, 0.0477063}; - static const float gaussrb[3] = {0.332406, 0.241376, 0.0924212}; + static const float gaussg[5] = {0.171582, 0.15839, 0.124594, 0.083518, 0.0477063};//sig=2.5 + static const float gaussrb[3] = {0.332406, 0.241376, 0.0924212};//sig=1.25 //char *buffer1; // vblsz*hblsz*3*2 //float (*blockshifts)[3][2]; // vblsz*hblsz*3*2 @@ -150,10 +150,10 @@ void RawImageSource::CA_correct_RT() { float (*rgb)[3]; // TS*TS*12 float (*grbdiff); // TS*TS*4 float (*gshift); // TS*TS*4 - float (*ghpfh); // TS*TS*4 - float (*ghpfv); // TS*TS*4 float (*rbhpfh); // TS*TS*4 float (*rbhpfv); // TS*TS*4 + float (*rblpfh); // TS*TS*4 + float (*rblpfv); // TS*TS*4 /* assign working space; this would not be necessary @@ -166,14 +166,16 @@ void RawImageSource::CA_correct_RT() { rgb = (float (*)[3]) buffer; grbdiff = (float (*)) (buffer + 12*TS*TS); gshift = (float (*)) (buffer + 16*TS*TS); - ghpfh = (float (*)) (buffer + 20*TS*TS); - ghpfv = (float (*)) (buffer + 24*TS*TS); - rbhpfh = (float (*)) (buffer + 28*TS*TS); - rbhpfv = (float (*)) (buffer + 32*TS*TS); + rbhpfh = (float (*)) (buffer + 20*TS*TS); + rbhpfv = (float (*)) (buffer + 24*TS*TS); + rblpfh = (float (*)) (buffer + 28*TS*TS); + rblpfv = (float (*)) (buffer + 32*TS*TS); - - vblsz=ceil((float)(height+border2)/(TS-border2))+2; - hblsz=ceil((float)(width+border2)/(TS-border2))+2; + if((height+border2)%(TS-border2)==0) vz1=1; else vz1=0; + if((width+border2)%(TS-border2)==0) hz1=1; else hz1=0; + + vblsz=ceil((float)(height+border2)/(TS-border2)+2+vz1); + hblsz=ceil((float)(width+border2)/(TS-border2)+2+hz1); /*buffer1 = (char *) malloc(4*vblsz*hblsz*3*2); merror(buffer1,"CA_correct()"); @@ -208,9 +210,8 @@ void RawImageSource::CA_correct_RT() { indx=row*width+col; indx1=rr*TS+cc; rgb[indx1][c] = (ri->data[row][col])/65535.0f; + //rgb[indx1][c] = image[indx][c]/65535.0f;//for dcraw implementation } - blockwt[vblock*hblsz+hblock]=0; - //blockwt[vblock*hblsz+hblock]=1; // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% //fill borders @@ -226,6 +227,7 @@ void RawImageSource::CA_correct_RT() { for (cc=ccmin; ccdata[(height-rr-2)][left+cc])/65535.0f; + //rgb[(rrmax+rr)*TS+cc][c] = (image[(height-rr-2)*width+left+cc][c])/65535.0f;//for dcraw implementation } } if (ccmin>0) { @@ -240,6 +242,7 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[(top+rr)][(width-cc-2)])/65535.0f; + //rgb[rr*TS+ccmax+cc][c] = (image[(top+rr)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation } } @@ -249,6 +252,7 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[border2-rr][border2-cc])/65535.0f; + //rgb[(rr)*TS+cc][c] = (rgb[(border2-rr)*TS+(border2-cc)][c]);//for dcraw implementation } } if (rrmaxdata[(height-rr-2)][(width-cc-2)])/65535.0f; + //rgb[(rrmax+rr)*TS+ccmax+cc][c] = (image[(height-rr-2)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation } } if (rrmin>0 && ccmaxdata[(border2-rr)][(width-cc-2)])/65535.0f; + //rgb[(rr)*TS+ccmax+cc][c] = (image[(border2-rr)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation } } if (rrmax0) { @@ -270,19 +276,22 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[(height-rr-2)][(border2-cc)])/65535.0f; + //rgb[(rrmax+rr)*TS+cc][c] = (image[(height-rr-2)*width+(border2-cc)][c])/65535.0f;//for dcraw implementation } } //end of border fill // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - for (j=0; j<2; j++) - for (k=0; k<3; k++) - for (c=0; c<3; c+=2) coeff[j][k][c]=0; + for (j=0; j<2; j++) + for (k=0; k<3; k++) + for (c=0; c<3; c+=2) { + coeff[j][k][c]=0; + } //end of initialization - + for (rr=3; rr < rr1-3; rr++) for (row=rr+top, cc=3, indx=rr*TS+cc; cc < cc1-3; cc++, indx++) { col = cc+left; @@ -294,7 +303,6 @@ void RawImageSource::CA_correct_RT() { wtd=1/SQR(eps+fabs(rgb[(rr-1)*TS+cc][1]-rgb[(rr+1)*TS+cc][1])+fabs(rgb[(rr)*TS+cc][c]-rgb[(rr+2)*TS+cc][c])+fabs(rgb[(rr+1)*TS+cc][1]-rgb[(rr+3)*TS+cc][1])); wtl=1/SQR(eps+fabs(rgb[(rr)*TS+cc+1][1]-rgb[(rr)*TS+cc-1][1])+fabs(rgb[(rr)*TS+cc][c]-rgb[(rr)*TS+cc-2][c])+fabs(rgb[(rr)*TS+cc-1][1]-rgb[(rr)*TS+cc-3][1])); wtr=1/SQR(eps+fabs(rgb[(rr)*TS+cc-1][1]-rgb[(rr)*TS+cc+1][1])+fabs(rgb[(rr)*TS+cc][c]-rgb[(rr)*TS+cc+2][c])+fabs(rgb[(rr)*TS+cc+1][1]-rgb[(rr)*TS+cc+3][1])); - //store in rgb array the interpolated G value at R/B grid points using directional weighted average rgb[indx][1]=(wtu*rgb[indx-v1][1]+wtd*rgb[indx+v1][1]+wtl*rgb[indx-1][1]+wtr*rgb[indx+1][1])/(wtu+wtd+wtl+wtr); @@ -302,46 +310,66 @@ void RawImageSource::CA_correct_RT() { if (row>-1 && row-1 && col0) { + CAshift[j][c]=coeff[j][1][c]/coeff[j][2][c]; + blockwt[vblock*hblsz+hblock]= areawt[j][c];//*coeff[j][2][c]/(eps+coeff[j][0][c]) ; + } else { + CAshift[j][c]=17.0; + blockwt[vblock*hblsz+hblock]=0; + } + + //CAshift[j][c]=coeff[j][1][c]/coeff[j][2][c]; + //blockwt[vblock*hblsz+hblock] = (float)(rr1-8)*(cc1-8)/4 * coeff[j][2][c]/(eps+coeff[j][0][c]) ; + //data structure = CAshift[vert/hor][color] //j=0=vert, 1=hor + if ((CAshift[j][c])<0) { GRBdir[j][c]=-1; } else { @@ -370,17 +407,15 @@ void RawImageSource::CA_correct_RT() { blocksqave[j][c] += SQR(CAshift[j][c]); blockdenom[j][c] += 1; } - }//vert/hor }//color - + /* CAshift[j][c] are the locations that minimize color difference variances; This is the approximate _optical_ location of the R/B pixels */ - - + for (c=0; c<3; c+=2) { //evaluate the shifts to the location that minimizes CA within the tile blockshifts[(vblock)*hblsz+hblock][c][0]=(CAshift[0][c]); //vert CA shift for R/B @@ -388,13 +423,19 @@ void RawImageSource::CA_correct_RT() { //data structure: blockshifts[blocknum][R/B][v/h] } } + //end of diagnostic pass for (j=0; j<2; j++) for (c=0; c<3; c+=2) { - blockvar[j][c] = blocksqave[j][c]/blockdenom[j][c]-SQR(blockave[j][c]/blockdenom[j][c]); + if (blockdenom[j][c]) { + blockvar[j][c] = blocksqave[j][c]/blockdenom[j][c]-SQR(blockave[j][c]/blockdenom[j][c]); + } else { + return; + } } - //end of diagnostic pass + //if (verbose) fprintf (stderr,_("tile variances %f %f %f %f \n"),blockvar[0][0],blockvar[1][0],blockvar[0][2],blockvar[1][2] ); + // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -447,17 +488,21 @@ void RawImageSource::CA_correct_RT() { if (p[4]<0) {GRBdir[dir][c]=-1;} else {GRBdir[dir][c]=1;} } + + //if (verbose) fprintf (stderr,_("tile vshift hshift (%d %d %4f %4f)...\n"),vblock, hblock, blockshifts[(vblock)*hblsz+hblock][c][0], blockshifts[(vblock)*hblsz+hblock][c][1]); + + //now prepare coefficient matrix - if (SQR(blockshifts[(vblock)*hblsz+hblock][c][0])>4*blockvar[0][c] || SQR(blockshifts[(vblock)*hblsz+hblock][c][1])>4*blockvar[1][c]) continue; + if (SQR(blockshifts[(vblock)*hblsz+hblock][c][0])>4.0*blockvar[0][c] || SQR(blockshifts[(vblock)*hblsz+hblock][c][1])>4.0*blockvar[1][c]) continue; numblox[c] += 1; for (dir=0; dir<2; dir++) { for (i=0; idata[row][col])/65535.0f; + //rgb[indx1][c] = image[indx][c]/65535.0f;//for dcraw implementation + if ((c&1)==0) rgb[indx1][1] = Gtmp[indx]; } // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -530,6 +586,8 @@ void RawImageSource::CA_correct_RT() { for (cc=ccmin; ccdata[(height-rr-2)][left+cc])/65535.0f; + //rgb[(rrmax+rr)*TS+cc][c] = (image[(height-rr-2)*width+left+cc][c])/65535.0f;//for dcraw implementation + rgb[(rrmax+rr)*TS+cc][1] = Gtmp[(height-rr-2)*width+left+cc]; } } @@ -546,6 +604,8 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[(top+rr)][(width-cc-2)])/65535.0f; + //rgb[rr*TS+ccmax+cc][c] = (image[(top+rr)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation + rgb[rr*TS+ccmax+cc][1] = Gtmp[(top+rr)*width+(width-cc-2)]; } } @@ -556,6 +616,8 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[border2-rr][border2-cc])/65535.0f; + //rgb[(rr)*TS+cc][c] = (rgb[(border2-rr)*TS+(border2-cc)][c]);//for dcraw implementation + rgb[(rr)*TS+cc][1] = Gtmp[(border2-rr)*width+border2-cc]; } } @@ -564,6 +626,8 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[(height-rr-2)][(width-cc-2)])/65535.0f; + //rgb[(rrmax+rr)*TS+ccmax+cc][c] = (image[(height-rr-2)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation + rgb[(rrmax+rr)*TS+ccmax+cc][1] = Gtmp[(height-rr-2)*width+(width-cc-2)]; } } @@ -572,6 +636,8 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[(border2-rr)][(width-cc-2)])/65535.0f; + //rgb[(rr)*TS+ccmax+cc][c] = (image[(border2-rr)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation + rgb[(rr)*TS+ccmax+cc][1] = Gtmp[(border2-rr)*width+(width-cc-2)]; } } @@ -580,6 +646,8 @@ void RawImageSource::CA_correct_RT() { for (cc=0; ccdata[(height-rr-2)][(border2-cc)])/65535.0f; + //rgb[(rrmax+rr)*TS+cc][c] = (image[(height-rr-2)*width+(border2-cc)][c])/65535.0f;//for dcraw implementation + rgb[(rrmax+rr)*TS+cc][1] = Gtmp[(height-rr-2)*width+(border2-cc)]; } } @@ -588,16 +656,6 @@ void RawImageSource::CA_correct_RT() { // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - blockshifts[(vblock)*hblsz+hblock][0][0] = blockshifts[(vblock)*hblsz+hblock][0][1] = 0; - blockshifts[(vblock)*hblsz+hblock][2][0] = blockshifts[(vblock)*hblsz+hblock][2][1] = 0; - for (i=0; ifabs(grbdiffint) ) { + rgb[indx][c]=RBint; + } } else { //gradient weights using difference from G at CA shift points and G at grid points - p[0]=1/(eps+fabs(rgb[(rr)*TS+cc][1]-gshift[(rr)*TS+cc])); - p[1]=1/(eps+fabs(rgb[(rr)*TS+cc][1]-gshift[(rr)*TS+cc-2*GRBdir[1][c]])); - p[2]=1/(eps+fabs(rgb[(rr)*TS+cc][1]-gshift[(rr-2*GRBdir[0][c])*TS+cc])); - p[3]=1/(eps+fabs(rgb[(rr)*TS+cc][1]-gshift[(rr-2*GRBdir[0][c])*TS+cc-2*GRBdir[1][c]])); + p[0]=1/(eps+fabs(rgb[indx][1]-gshift[indx])); + p[1]=1/(eps+fabs(rgb[indx][1]-gshift[indx-2*GRBdir[1][c]])); + p[2]=1/(eps+fabs(rgb[indx][1]-gshift[(rr-2*GRBdir[0][c])*TS+cc])); + p[3]=1/(eps+fabs(rgb[indx][1]-gshift[(rr-2*GRBdir[0][c])*TS+cc-2*GRBdir[1][c]])); - grbdiffint = (p[0]*grbdiff[(rr)*TS+cc]+p[1]*grbdiff[(rr)*TS+cc-2*GRBdir[1][c]]+ \ + grbdiffint = (p[0]*grbdiff[indx]+p[1]*grbdiff[indx-2*GRBdir[1][c]]+ \ p[2]*grbdiff[(rr-2*GRBdir[0][c])*TS+cc]+p[3]*grbdiff[(rr-2*GRBdir[0][c])*TS+cc-2*GRBdir[1][c]])/(p[0]+p[1]+p[2]+p[3]); - //there is an assumption that the grid point is no more than 2 pixels from the optical point; perhaps should correct for this??? //now determine R/B at grid points using interpolated color differences and interpolated G value at grid point - grbdiffold = rgb[(rr)*TS+cc][1]-rgb[(rr)*TS+cc][c]; - if (fabs(grbdiffold)>fabs(grbdiffint) ) { - rgb[(rr)*TS+cc][c]=rgb[(rr)*TS+cc][1]-grbdiffint; + rgb[indx][c]=rgb[indx][1]-grbdiffint; } - if (grbdiffold*grbdiffint<0) { - rgb[(rr)*TS+cc][c]=rgb[(rr)*TS+cc][1]; - } - + } + + //if color difference interpolation overshot the correction, just desaturate + if (grbdiffold*grbdiffint<0) { + rgb[indx][c]=rgb[indx][1]-0.5*(grbdiffold+grbdiffint); } } @@ -673,8 +733,11 @@ void RawImageSource::CA_correct_RT() { c = FC(row,col); ri->data[row][col] = CLIP((int)(65535.0f*rgb[(rr)*TS+cc][c] + 0.5f)); + //image[indx][c] = CLIP((int)(65535.0*rgb[(rr)*TS+cc][c] + 0.5));//for dcraw implementation } + + if(plistener) plistener->setProgress(fabs((float)top/height)); } // clean up @@ -685,9 +748,8 @@ void RawImageSource::CA_correct_RT() { #undef TS -#undef polyord -#undef numpar -#undef border +//#undef border +//#undef border2 #undef PIX_SORT #undef SQR diff --git a/rtengine/amaze_interpolate_RT.cc b/rtengine/amaze_interpolate_RT.cc index 64f9d0e41..47770a3c1 100644 --- a/rtengine/amaze_interpolate_RT.cc +++ b/rtengine/amaze_interpolate_RT.cc @@ -64,6 +64,7 @@ void RawImageSource::amaze_demosaic_RT() { static const float arthresh=0.75; static const float nyqthresh=0.5;//0.5 static const float pmthresh=0.25;//0.25 + static const float lbd=0.66, ubd=1.5; static const float gaussodd[4] = {0.14659727707323927f, 0.103592713382435f, 0.0732036125103057f, 0.0365543548389495f};//gaussian on 5x5 quincunx, sigma=1.2 float gaussgrad[6] = {0.07384411893421103f, 0.06207511968171489f, 0.0521818194747806f, \ @@ -78,7 +79,9 @@ void RawImageSource::amaze_demosaic_RT() { int bottom, right, row, col; int rr, cc, rr1, cc1, c, indx, indx1, dir, i, j, sgn; - + float crse, crnw, crne, crsw, rbse, rbnw, rbne, rbsw, wtse, wtnw, wtsw, wtne; + float pmwtalt; + float cru, crd, crl, crr; float vwt, hwt, Gintv, Ginth; float guar, gdar, glar, grar, guha, gdha, glha, grha, Ginthar, Ginthha, Gintvar, Gintvha, hcdaltvar, vcdaltvar; @@ -86,7 +89,7 @@ void RawImageSource::amaze_demosaic_RT() { float sumh, sumv, sumsqh, sumsqv, areawt; float nyqtest, vcdvar, hcdvar, hvwtalt, vo, ve, gradp, gradm, gradv, gradh, gradpm, gradhv; float vcdvar1, hcdvar1, varwt, diffwt; - float rbvarp, rbvarm, crp, crm, rbp, rbm; + float rbvarp, rbvarm; float gu, gd, gl, gr; float gvarh, gvarv; float g[4], f[4]; @@ -124,14 +127,17 @@ void RawImageSource::amaze_demosaic_RT() { float (*Dgrbpsq1); // TS*TS*4 float (*Dgrbmsq1); // TS*TS*4 float (*cfa); // TS*TS*4 + float (*pmwt); // TS*TS*4 + float (*rbp); // TS*TS*4 + float (*rbm); // TS*TS*4 - int (*nyquist); // TS*TS*4 + int (*nyquist); // TS*TS*4 // assign working space - buffer = (char *) malloc(144*TS*TS); + buffer = (char *) malloc(156*TS*TS); //merror(buffer,"amaze_interpolate()"); - memset(buffer,0,144*TS*TS); + memset(buffer,0,156*TS*TS); // rgb array rgb = (float (*)[3]) buffer; //pointers to array delh = (float (*)) (buffer + 12*TS*TS); @@ -164,8 +170,11 @@ void RawImageSource::amaze_demosaic_RT() { Dgrbpsq1 = (float (*)) (buffer + 128*TS*TS); Dgrbmsq1 = (float (*)) (buffer + 132*TS*TS); cfa = (float (*)) (buffer + 136*TS*TS); + pmwt = (float (*)) (buffer + 140*TS*TS); + rbp = (float (*)) (buffer + 144*TS*TS); + rbm = (float (*)) (buffer + 148*TS*TS); - nyquist = (int (*)) (buffer + 140*TS*TS); + nyquist = (int (*)) (buffer + 152*TS*TS); // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -227,6 +236,8 @@ void RawImageSource::amaze_demosaic_RT() { indx=row*width+col; indx1=rr*TS+cc; rgb[indx1][c] = (ri->data[row][col])/65535.0f; + //rgb[indx1][c] = image[indx][c]/65535.0f;//for dcraw implementation + cfa[indx1] = rgb[indx1][c]; } // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -236,7 +247,6 @@ void RawImageSource::amaze_demosaic_RT() { for (cc=ccmin; ccdata[(height-rr-2)][left+cc])/65535.0f; - cfa[(rrmax+rr)*TS+cc] = rgb[(rrmax+rr)*TS+cc][c]; + //rgb[(rrmax+rr)*TS+cc][c] = (image[(height-rr-2)*width+left+cc][c])/65535.0f;//for dcraw implementation } } if (ccmin>0) { @@ -252,7 +262,6 @@ void RawImageSource::amaze_demosaic_RT() { for (cc=0; cc<16; cc++) { c=FC(rr,cc); rgb[rr*TS+cc][c] = rgb[rr*TS+32-cc][c]; - cfa[rr*TS+cc] = rgb[rr*TS+cc][c]; } } if (ccmaxdata[(top+rr)][(width-cc-2)])/65535.0f; - cfa[rr*TS+ccmax+cc] = rgb[rr*TS+ccmax+cc][c]; + //rgb[rr*TS+ccmax+cc][c] = (image[(top+rr)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation } } - //also fill the image corners + //also, fill the image corners if (rrmin>0 && ccmin>0) { for (rr=0; rr<16; rr++) for (cc=0; cc<16; cc++) { c=FC(rr,cc); rgb[(rr)*TS+cc][c] = (ri->data[32-rr][32-cc])/65535.0f; - cfa[(rr)*TS+cc] = rgb[(rr)*TS+cc][c]; + //rgb[(rr)*TS+cc][c] = (rgb[(32-rr)*TS+(32-cc)][c]);//for dcraw implementation } } if (rrmaxdata[(height-rr-2)][(width-cc-2)])/65535.0f; - cfa[(rrmax+rr)*TS+ccmax+cc] = rgb[(rrmax+rr)*TS+ccmax+cc][c]; + //rgb[(rrmax+rr)*TS+ccmax+cc][c] = (image[(height-rr-2)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation } } if (rrmin>0 && ccmaxdata[(rr)][(width-cc-2)])/65535.0f; - cfa[(rr)*TS+ccmax+cc] = rgb[(rr)*TS+ccmax+cc][c]; + rgb[(rr)*TS+ccmax+cc][c] = (ri->data[(32-rr)][(width-cc-2)])/65535.0f; + //rgb[(rr)*TS+ccmax+cc][c] = (image[(32-rr)*width+(width-cc-2)][c])/65535.0f;//for dcraw implementation } } if (rrmax0) { for (rr=0; rr<16; rr++) for (cc=0; cc<16; cc++) { c=FC(rr,cc); - rgb[(rrmax+rr)*TS+cc][c] = (ri->data[(height-rr-2)][cc])/65535.0f; - cfa[(rrmax+rr)*TS+cc] = rgb[(rrmax+rr)*TS+cc][c]; + rgb[(rrmax+rr)*TS+cc][c] = (ri->data[(height-rr-2)][(32-cc)])/65535.0f; + //rgb[(rrmax+rr)*TS+cc][c] = (image[(height-rr-2)*width+(32-cc)][c])/65535.0f;//for dcraw implementation } } - + //end of border fill // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -349,10 +358,10 @@ void RawImageSource::amaze_demosaic_RT() { rbint[indx]=0; //color ratios in each cardinal direction - cru = cfa[indx-v1]*(eps+dirwts[indx-v2][0]+dirwts[indx][0])/(eps+dirwts[indx-v2][0]*cfa[indx]+dirwts[indx][0]*cfa[indx-v2]); - crd = cfa[indx+v1]*(eps+dirwts[indx+v2][0]+dirwts[indx][0])/(eps+dirwts[indx+v2][0]*cfa[indx]+dirwts[indx][0]*cfa[indx+v2]); - crl = cfa[indx-1]*(eps+dirwts[indx-2][1]+dirwts[indx][1])/(eps+dirwts[indx-2][1]*cfa[indx]+dirwts[indx][1]*cfa[indx-2]); - crr = cfa[indx+1]*(eps+dirwts[indx+2][1]+dirwts[indx][1])/(eps+dirwts[indx+2][1]*cfa[indx]+dirwts[indx][1]*cfa[indx+2]); + cru = cfa[indx-v1]*(dirwts[indx-v2][0]+dirwts[indx][0])/(dirwts[indx-v2][0]*cfa[indx]+dirwts[indx][0]*cfa[indx-v2]); + crd = cfa[indx+v1]*(dirwts[indx+v2][0]+dirwts[indx][0])/(dirwts[indx+v2][0]*cfa[indx]+dirwts[indx][0]*cfa[indx+v2]); + crl = cfa[indx-1]*(dirwts[indx-2][1]+dirwts[indx][1])/(dirwts[indx-2][1]*cfa[indx]+dirwts[indx][1]*cfa[indx-2]); + crr = cfa[indx+1]*(dirwts[indx+2][1]+dirwts[indx][1])/(dirwts[indx+2][1]*cfa[indx]+dirwts[indx][1]*cfa[indx+2]); guha=cfa[indx-v1]+0.5*(cfa[indx]-cfa[indx-v2]); gdha=cfa[indx+v1]+0.5*(cfa[indx]-cfa[indx+v2]); @@ -398,17 +407,21 @@ void RawImageSource::amaze_demosaic_RT() { if (hcdaltvar3) pmwt[indx]=1; + //or not + if (areawt<2) pmwt[indx]=0; + }*/ + for (rr=8; rr (0.33*(2*Ginth+rbint[indx]))) Ginth=ULIM(Ginth,cfa[indx-1],cfa[indx+1]); - if (rbint[indx]-2*Gintv > (0.33*(2*Gintv+rbint[indx]))) Gintv=ULIM(Gintv,cfa[indx-v1],cfa[indx+v1]); + //if (rbint[indx]-2*Ginth > (0.33*(2*Ginth+rbint[indx]))) Ginth=ULIM(Ginth,cfa[indx-1],cfa[indx+1]); + //if (rbint[indx]-2*Gintv > (0.33*(2*Gintv+rbint[indx]))) Gintv=ULIM(Gintv,cfa[indx-v1],cfa[indx+v1]); + Ginth=LIM(Ginth, lbd*MIN(cfa[indx-1],cfa[indx+1]), ubd*MAX(cfa[indx-1],cfa[indx+1])); + Gintv=LIM(Gintv, lbd*MIN(cfa[indx-v1],cfa[indx+v1]), ubd*MAX(cfa[indx-v1],cfa[indx+v1])); + //Galt[indx] = Ginth*(1-hvwt[indx]) + Gintv*hvwt[indx]; rgb[indx][1] = Ginth*(1-hvwt[indx]) + Gintv*hvwt[indx]; + Dgrb[indx][0] = rgb[indx][1]-cfa[indx]; + + rgb[indx][2-FC(rr,cc)]=2*rbint[indx]-cfa[indx]; } //end of diagonal interpolation correction //t2_diag += clock() - t1_diag; @@ -744,6 +829,12 @@ void RawImageSource::amaze_demosaic_RT() { green[row][col] = CLIP((int)(65535.0f*rgb[indx][1] + 0.5f)); blue[row][col] = CLIP((int)(65535.0f*rgb[indx][2] + 0.5f)); + //for dcraw implementation + //for (c=0; c<3; c++){ + // image[indx][c] = CLIP((int)(65535.0f*rgb[rr*TS+cc][c] + 0.5f)); + //} + + } //end of main loop