modifying convolution in LMA (24c49504) · Commits · Elphel / imagej-elphel

src/main/java/EyesisDCT.java

+40 −54

Original line number	Diff line number	Diff line
		@@ -214,7 +214,7 @@ public class EyesisDCT {
		for (int i = 0; i <target_kernel.length; i++){
		target_kernel[i]/=s;
		}
		if (globalDebugLevel>0){
		if (globalDebugLevel > 1){ // was already close to 1.0
		System.out.println(tileX+"/"+tileY+ " s="+s);
		}
		}
		@@ -634,53 +634,6 @@ public class EyesisDCT {
		for (int nc = 0; nc < nColors; nc++){
		for (int tileY = 0; tileY < numVert; tileY++){
		for (int tileX = 0; tileX < numHor; tileX++){
		// extract DCT (symmetrical) kernels
		int sym_kernel_start_index = (sym_kernel_inc_index * tileY + tileX) * dct_size;
		for (int i = 0; i < dct_size;i++){
		System.arraycopy( // copy one kernel line
		kernels[chn].sym_kernels[nc],
		sym_kernel_start_index + i * sym_kernel_inc_index,
		kernels[chn].st_kernels[nc][tileY][tileX],
		i * dct_size,
		dct_size);
		}
		double scale = dct_size;
		if (dct_parameters.normalize) {
		double sum = 0.0;
		for (int i=0; i < dct_size;i++) {
		for (int j=0; j < dct_size;j++) {
		double d = kernels[chn].st_kernels[nc][tileY][tileX][i*dct_size+j];
		if (i > 0) d *= 2.0;
		if (j > 0) d *= 2.0;
		sum += d;
		}
		}

		scale /= sum;
		if ((tileY == ((dct_parameters.tileY+1)/2)) && (tileX ==((dct_parameters.tileX+1)/2)) && (nc == 2)) {
		System.out.println("(kernel) tileY="+tileY+"(kernel) tileX="+tileX+" sum="+sum + " scale="+scale);
		}

		}
		for (int i=0; i < kernels[chn].st_kernels[nc][tileY][tileX].length;i++) {
		kernels[chn].st_kernels[nc][tileY][tileX][i] *= scale;
		}
		// Make a copy of direct kernels (debug feature, may be removed later)
		for (int i = 0; i < dct_size;i++){
		System.arraycopy( // copy one kernel line
		kernels[chn].st_kernels[nc][tileY][tileX],
		i * dct_size,
		kernels[chn].st_direct[nc][tileY][tileX],
		i * dct_size,
		dct_size);
		}
		// scale st_direct back to ~original
		for (int i=0; i < kernels[chn].st_kernels[nc][tileY][tileX].length;i++) {
		kernels[chn].st_direct[nc][tileY][tileX][i] /= dct_size;
		}

		kernels[chn].st_kernels[nc][tileY][tileX]= dtt.dttt_iii(kernels[chn].st_kernels[nc][tileY][tileX]);

		// extract asymmetrical kernel and convert it to list of values and indices (as arrays as the length is known)
		int asym_kernel_start_index = (asym_kernel_inc_index * tileY + tileX)* asym_size;
		int indx = 0;
		@@ -705,26 +658,59 @@ public class EyesisDCT {
		}
		}

		double scale_asym = 0.0;
		if (dct_parameters.normalize) {
		double sum = 0.0;
		for (int i = 0; i < kernels[chn].asym_val[nc][tileY][tileX].length;i++){
		sum += kernels[chn].asym_val[nc][tileY][tileX][i];
		scale_asym += kernels[chn].asym_val[nc][tileY][tileX][i];
		if ((tileY==67) && (tileX==125)) {
		System.out.println("i="+i+", sum="+sum);
		System.out.println("i="+i+", sum="+scale_asym);
		}
		}
		for (int i = 0; i < kernels[chn].asym_val[nc][tileY][tileX].length;i++){
		kernels[chn].asym_val[nc][tileY][tileX][i] /= sum;
		kernels[chn].asym_val[nc][tileY][tileX][i] /= scale_asym;
		}
		if ((tileY==67) && (tileX==125)) {
		System.out.println("sum="+sum+", normalized:");
		System.out.println("sum="+scale_asym+", normalized:");

		for (int i=0; i<kernels[chn].asym_indx[nc][tileY][tileX].length; i++){
		System.out.println("kernels["+chn+"].asym_val["+nc+"]["+tileY+"]["+tileX+"]["+i+"]="+kernels[chn].asym_val[nc][tileY][tileX][i]);
		System.out.println("kernels["+chn+"].asym_indx["+nc+"]["+tileY+"]["+tileX+"]["+i+"]="+kernels[chn].asym_indx[nc][tileY][tileX][i]);
		}
		}
		} else {
		scale_asym = 1.0;
		}

		// extract DCT (symmetrical) kernels
		int sym_kernel_start_index = (sym_kernel_inc_index * tileY + tileX) * dct_size;
		for (int i = 0; i < dct_size;i++){
		System.arraycopy( // copy one kernel line
		kernels[chn].sym_kernels[nc],
		sym_kernel_start_index + i * sym_kernel_inc_index,
		kernels[chn].st_kernels[nc][tileY][tileX],
		i * dct_size,
		dct_size);
		}
		if (scale_asym != 1.0){
		for (int i=0; i < kernels[chn].st_kernels[nc][tileY][tileX].length;i++) {
		kernels[chn].st_kernels[nc][tileY][tileX][i] *= scale_asym;
		}
		}
		// Make a copy of direct kernels (debug feature, may be removed later)
		for (int i = 0; i < dct_size;i++){
		System.arraycopy( // copy one kernel line
		kernels[chn].st_kernels[nc][tileY][tileX],
		i * dct_size,
		kernels[chn].st_direct[nc][tileY][tileX],
		i * dct_size,
		dct_size);
		}
		// scale so multiplication will not change normalization
		for (int i=0; i < kernels[chn].st_kernels[nc][tileY][tileX].length;i++) {
		kernels[chn].st_kernels[nc][tileY][tileX][i] *= dct_size;
		}

		kernels[chn].st_kernels[nc][tileY][tileX]= dtt.dttt_iii(kernels[chn].st_kernels[nc][tileY][tileX]);
		}
		// System.out.println("tileY="+tileY);
		}

src/main/java/FactorConvKernel.java

+37 −16

Original line number	Diff line number	Diff line
		@@ -487,6 +487,7 @@ public class FactorConvKernel {
		int conv_size = asym_size + 2*sym_radius-2;
		int conv_len = conv_size * conv_size;
		int sym_rad_m1 = sym_radius - 1; // 7
		int sym_rad2 = 2*sym_radius; // 16
		double [] fX = new double [justConvolved? conv_len: this.weight.length];
		// calculate convolution, for kernels - regardless of kernels enabled/disabled
		// calculate convolution part
		@@ -500,15 +501,25 @@ public class FactorConvKernel {
		for (int ai = 0; ai < asym_size; ai ++){
		int si = (ci - ai) - sym_rad_m1;
		if (si < 0) si = -si;
		if (si < sym_radius) {
		int sgni = 1;
		if (si > sym_radius) {
		sgni = -1;
		si = sym_rad2 - si;
		}
		if (si < sym_radius) { // skip si == sym_radius, coefficient is 0 (WA)
		for (int aj = 0; aj < asym_size; aj ++){
		int aindx = ai*asym_size + aj;
		if (!skip_disabled_asym \|\| kernel_masks[1][aindx]){
		int sj = (cj - aj) - sym_rad_m1;
		if (sj < 0) sj = -sj;
		if (sj < sym_radius) {
		int sgn = sgni;
		if (sj > sym_radius) {
		sgn = -sgn;
		sj = sym_rad2 - sj;
		}
		if (sj < sym_radius) { // skip sj == sym_radius, coefficient is 0 (WA)
		int sindx = si * sym_radius + sj;
		fX[cindx] += kernels[0][sindx] * kernels[1][aindx];
		fX[cindx] += sgnkernels[0][sindx] kernels[1][aindx];
		}
		}
		}
		@@ -551,6 +562,8 @@ public class FactorConvKernel {
		int conv_size = asym_size + 2*sym_radius-2;
		int conv_len = conv_size * conv_size;
		int sym_rad_m1 = sym_radius - 1; // 7
		int sym_rad2 = 2*sym_radius; // 16
		int sym_rad4 = 4*sym_radius; // 32
		// calculate convolution part
		for (int ci =0; ci < conv_size; ci++) for (int cj =0; cj < conv_size; cj++){
		int cindx = ci*conv_size + cj;
		@@ -559,20 +572,37 @@ public class FactorConvKernel {
		for (int ai = 0; ai < asym_size; ai ++){
		int si = (ci - ai) - sym_rad_m1;
		if (si < 0) si = -si;
		if (si < sym_radius) {
		int sgni = 1;
		if (si > sym_rad2) si = sym_rad4 - si;
		if (si > sym_radius) {
		sgni = -1;
		si = sym_rad2 - si;
		}
		if (si < sym_radius) { // skip si == sym_radius, coefficient is 0 (WA)
		for (int aj = 0; aj < asym_size; aj ++){
		int aindx = ai*asym_size + aj;
		int apar_indx = map_to_pars[1][aindx];
		int sj = (cj - aj) - sym_rad_m1;
		if (sj < 0) sj = -sj;
		if (sj < sym_radius) {
		int sgn = sgni;
		if (sj > sym_rad2) sj = sym_rad4 - sj;
		if (sj > sym_radius) {
		sgn = -sgn;
		sj = sym_rad2 - sj;
		}
		if (sj < sym_radius) { // skip sj == sym_radius, coefficient is 0 (WA)
		int sindx = si * sym_radius + sj;
		// if (sindx <0){
		// System.out.println(
		// "ci="+ci+" cj="+cj+" si="+si+" sj="+sj+" ai="+ai+" aj="+aj+
		// " sgni="+sgni+" sgn="+sgn+" sym_rad2="+sym_rad2);
		// }
		if (apar_indx >= 0){
		jacobian[apar_indx][cindx] += kernels[0][sindx];
		jacobian[apar_indx][cindx] += sgn*kernels[0][sindx];
		}
		int spar_indx = map_to_pars[0][sindx];
		if ((spar_indx>=0) && (!skip_disabled_asym \|\| kernel_masks[1][aindx])){
		jacobian[spar_indx][cindx] += kernels[1][aindx];
		jacobian[spar_indx][cindx] += sgn*kernels[1][aindx];
		}
		}
		}
		@@ -937,17 +967,8 @@ public class FactorConvKernel {
		this.target_kernel = target_kernel;
		this.asym_pixels = asym_pixels;
		this.asym_distance = asym_distance;
		/*
		double s = 0.0;
		for (int i = 0; i<target_kernel.length; i++){
		s+= target_kernel[i]*target_kernel[i];
		}
		this.goal_rms_pure = Math.sqrt(s/target_kernel.length)*fact_precision;
		*/

		double [] RMSes = null;
		// double [] bestRms = new double [asym_pixels];
		// int [] enPixels = new int [asym_pixels];
		this.startTime=System.nanoTime(); // need local?
		int numWeakest = 0;
		int numAny=0;

src/main/java/ImageDtt.java

+23 −0

Original line number	Diff line number	Diff line
		@@ -350,6 +350,8 @@ public class ImageDtt {
		int ady = (dy>=0)?dy:-dy;
		int sgny = 1;
		int y = i - dy;
		/*

		if (y < 0){
		y = -1 -y;
		sgny = -sgny;
		@@ -366,10 +368,21 @@ public class ImageDtt {
		sgny = 0;
		}
		}

		*/
		if (y < 0){
		y = -1 -y;
		sgny = -sgny;
		}
		if (y >= n2){
		y = 2*n2 - y -1;
		sgny = -sgny;
		}
		for (int dx = -dct_size +1; dx < dct_size; dx++){
		int adx = (dx >= 0)? dx:-dx;
		int sgn = sgny;
		int x = j - dx;
		/*
		if (x < 0){
		x = -1 -x;
		sgn = -sgn;
		@@ -386,6 +399,16 @@ public class ImageDtt {
		sgn = 0;
		}
		}

		*/
		if (x < 0){
		x = -1 -x;
		sgn = -sgn;
		}
		if (x >= n2){
		x = 2*n2 - x -1;
		sgn = -sgn;
		}
		sym_conv[indx] += sgndir_sym[ady dct_size + adx] * tile_in[y * n2 + x];
		s0+=dir_sym[ady * dct_size + adx];
		if ((tileY == debug_tileY) && (tileX == debug_tileX) && (color == 2) &&

Admin message