Commit e6a1d91b authored by Andrey Filippov's avatar Andrey Filippov
Browse files

center-of-mass for 16, propagated changes (4->16) for disparity_map

parent d8b8ffdc

src/main/java/com/elphel/imagej/tileprocessor/Correlation2d.java

+94 −59
Original line number Diff line number Diff line
@@ -1971,6 +1971,7 @@ public class Correlation2d { 
		}

		return getMaxXYInt( // find integer pair or null if below threshold // USED in lwir

				data,      // [data_size * data_size]

				null,

				data_width,

				center_row,

				axis_only,
@@ -1992,12 +1993,14 @@ public class Correlation2d { 
	

	public int [] getMaxXYInt( // find integer pair or null if below threshold // USED in lwir

			double [] data,      // [data_size * data_size]

			double [] disp_str,  // if not null, will return {disparity, strength}

			int       data_width,

			int       center_row,

			boolean   axis_only,

			double    minMax,    // minimal value to consider (at integer location, not interpolated)

			boolean   debug)

	{

		//mcorr_comb_disp

//		int data_width = 2 * transform_size - 1;

//		int data_height = data.length / data_width;

		int center = data_width / 2; // transform_size - 1;
@@ -2034,6 +2037,10 @@ public class Correlation2d { 
		if (debug){

			System.out.println("getMaxXYInt() -> "+rslt[0]+"/"+rslt[1]);

		}

		if (disp_str != null) {

			disp_str[0] = -rslt[0] * mcorr_comb_disp;

			disp_str[1] = data[imx];

		}

		return rslt;

	}
@@ -2108,10 +2115,14 @@ public class Correlation2d { 
	 */

	public double [] getMaxXCm( // get fractional center as a "center of mass" inside circle/square from the integer max // USED in lwir

			double [] data,      // [data_size * data_size]

			int       data_width,      //  = 2 * transform_size - 1;

			int       center_row,

			int       ixcenter,  // integer center x

			boolean   debug) {

		return getMaxXCm( // get fractional center as a "center of mass" inside circle/square from the integer max

				data,      // double [] data,      // [data_size * data_size]

				data_width,           // int       data_width,      //  = 2 * transform_size - 1;

				center_row,           // int       center_row,

				ixcenter,  // int       ixcenter,  // integer center x

				this.corr_wndy, // double [] window_y,  // (half) window function in y-direction(perpendicular to disparity: for row0  ==1

				this.corr_wndx, // double [] window_x,  // half of a window function in x (disparity) direction
@@ -2120,27 +2131,63 @@ public class Correlation2d { 
	

	public double [] getMaxXCmNotch( // get fractional center as a "center of mass" inside circle/square from the integer max // not used in lwir

			double [] data,      // [data_size * data_size]

			int       data_width,      //  = 2 * transform_size - 1;

			int       center_row,

			int       ixcenter,  // integer center x

			boolean   debug) {

		return getMaxXCm(             // get fractional center as a "center of mass" inside circle/square from the integer max

				data,                 // double [] data,      // [data_size * data_size]

				data_width,           // int       data_width,      //  = 2 * transform_size - 1;

				center_row,           // int       center_row,

				ixcenter,             // int       ixcenter,  // integer center x

				this.corr_wndy_notch, // double [] window_y,  // (half) window function in y-direction(perpendicular to disparity: for row0  ==1

				this.corr_wndx,       // double [] window_x,  // half of a window function in x (disparity) direction

				debug);               // boolean   debug);

	}

	

	@Deprecated

	public double [] getMaxXCm( // get fractional center as a "center of mass" inside circle/square from the integer max // not used in lwir

			double [] data,      // [data_size * data_size]

			int       ixcenter,  // integer center x

			boolean   debug) {

		return getMaxXCm( // get fractional center as a "center of mass" inside circle/square from the integer max // not used in lwir

				data,      // [data_size * data_size]

				2 * transform_size - 1, // int       data_width, 

				transform_size - 1,     // int       center_row,

				ixcenter,  // integer center x

				debug);

	}

	

	@Deprecated

	public double [] getMaxXCmNotch( // get fractional center as a "center of mass" inside circle/square from the integer max // not used in lwir

			double [] data,      // [data_size * data_size]

			int       ixcenter,  // integer center x

			boolean   debug) {

		return getMaxXCmNotch( // get fractional center as a "center of mass" inside circle/square from the integer max // not used in lwir

				data,      // [data_size * data_size]

				2 * transform_size - 1, // int       data_width, 

				transform_size - 1,     // int       center_row,

				ixcenter,  // integer center x

				debug);

	}

	

	// No shift by 0.5 for 2021

	public double [] getMaxXCm( // get fractional center as a "center of mass" inside circle/square from the integer max // USED in lwir

			double [] data,      // rectangular strip of 1/2 of the correlation are with odd rows shifted by 1/2 pixels

			int       data_width,      //  = 2 * transform_size - 1;

			int       center_row,

			int       ixcenter,  // integer center x

			double [] window_y,  // (half) window function in y-direction(perpendicular to disparity: for row0  ==1

			double [] window_x,  // half of a window function in x (disparity) direction

			boolean   debug) {

		int center = transform_size - 1;

		int data_width = 2 * transform_size - 1;

		int data_height = data.length/data_width;

		double wy_scale = 1.0;

		int center_x = (data_width - 1)/2; //  = transform_size - 1;

		int x0 = center_x + ixcenter; // index of the argmax, starting with 0



//		int data_width = 2 * transform_size - 1;

//		int data_height = data.length/data_width;

//		double wy_scale = 1.0;

		/*

		if (data_height > window_y.length) {

			data_height = window_y.length;

		} else if (data_height < window_y.length) { // re-  // not used in lwir
@@ -2148,70 +2195,58 @@ public class Correlation2d { 
			for (int i = 1; i < data_height; i++) swy += window_y[i];

			wy_scale = 1.0/swy;

		}



		double [][]dbg_data = null;

		if (debug) {

			String [] dbg_titles = {"strip","*wnd_y"};

			dbg_data = new double [2][];

			dbg_data[0] =  debugStrip3(data);

			double [] data_0 = data.clone();

	    */

		double w_scale = 1.0;

		if ((center_row + window_y.length > data_height) || (center_row - window_y.length < 0)) {

			double sw = 0.0;

			for (int i = 0; i < data_height; i++) {

				for (int j = 0; j <  data_width; j++) {

					data_0[i * data_width + j] *= (i < window_y.length) ? (wy_scale * window_y[i]): 0.0;

				int dy = i - center_row;

				int ady = (dy > 0) ? dy : -dy;

				if (ady < window_y.length) {

					sw += window_y[ady];

				}

			}

			dbg_data[1] =  debugStrip3(data_0);

			int long_width = 2 * (2 * transform_size-1);

			if (dbg_data[0] != null) {

				(new ShowDoubleFloatArrays()).showArrays(

						dbg_data,

						long_width,

						dbg_data[0].length/long_width,

						true,

						"Strip",

						dbg_titles);

			w_scale /= sw; 

		}



			System.out.println("getMaxXCm(), ixcenter = "+ixcenter);

			for (int dy = 0; dy < data_height; dy++) {

				if ((dy & 1) != 0) System.out.print("    ");

				for (int dx = 0; dx < data_width; dx++) {

					System.out.print(String.format(" %8.5f", data[dy * data_width + dx]));

		if ((x0 + window_x.length > data_width) || (x0 - window_x.length < 0)) {

			double sw = 0.0;

			for (int i = 0; i < data_width; i++) {

				int dx = i - x0;

				int adx = (dx > 0) ? dx : -dx;

				if (adx < window_x.length) {

					sw += window_x[adx];

				}

				System.out.println();

			}

			System.out.println();

			w_scale /= sw; 

		}

		

		

		//		double [][]dbg_data = null;

		double s0=0.0, sx=0.0, sx2 = 0.0;

		int x0 = center + ixcenter; // index of the argmax, starting with 0

		for (int dy = 0; dy < data_height; dy++) {

			int odd = dy & 1;

			double wy = ((dy == 0)? wy_scale: (2.0 * wy_scale))*window_y[dy];

			int indx0 = data_width * dy;

			for (int adx = odd; adx < window_x.length; adx+=2) { // index in window_x

				for (int dir = (adx == 0)?1:-1; dir <= 1; dir+=2) {

					// calculate data index

					int idx = (adx * dir) >> 1;

					int x = 2 * idx + odd;

					int x1 = x0 + idx; // correct

					if (debug)	System.out.print(String.format(" %2d:%2d:%d %3d", dy,adx,dir,x));

					if ((x1 >= 0 ) && (x1 < data_width)) {

						double d = data[indx0+x1];

///						if (!Double.isNaN(d)) {

		for (int iy = 0; iy < data_height; iy++) {

			int dy = iy - center_row;

			int ady = (dy > 0) ? dy : -dy;

			if (ady < window_y.length) {

				double wy = w_scale * window_y[ady];

				int indx0 = data_width * iy;

				for (int ix = 0; ix < data_width; ix++) {

					int dx = ix - x0; // 0 at argmax

					int adx = (dx > 0) ? dx : -dx;

					if (adx < window_x.length) {

						double d = data[indx0 + ix];

						if (debug)	System.out.print(String.format(" %2d:%2d:%8.5f ", dy, dx, d));

						if (!Double.isNaN(d) && (d > 0.0)) { // with negative d s0 can get very low value (or even negative)

							d*= wy*window_x[adx];

							double w = wy * window_x[adx];

							d*=    w;

							s0+=   d;

							sx +=  d * x; // result x is twice larger (corresponds to window_x)

							sx2 += d * x * x;

							if (debug)	System.out.print(String.format("%8.5f", data[indx0+x1])); //d));

							sx +=  d * dx;

							sx2 += d * dx * dx;

						}						

					} else {

						if (debug)	System.out.print("********");

					}

				}

			}

			if (debug)	System.out.println();



		}

		if (debug){

			System.out.println("getMaxXCm() -> s0="+s0+", sx="+sx+", sx2="+sx2+", ixcenter="+ixcenter);
@@ -2220,16 +2255,16 @@ public class Correlation2d { 
		if (s0 == 0.0) return null;



		double [] rslt = {

				ixcenter + sx/s0/2,              // new center in disparity units, relative to the correlation center

				(ixcenter + sx/s0)* mcorr_comb_disp, // /2,              // new center in disparity units, relative to the correlation center

				s0,                                  // total "weight"

				Math.sqrt(s0*sx2 - sx*sx)/s0/2}; // standard deviation in disparity units (divide weight by the standard deviation for quality?)

				(Math.sqrt(s0*sx2 - sx*sx)/s0)* mcorr_comb_disp}; // /2}; // standard deviation in disparity units (divide weight by the standard deviation for quality?)

		if (debug){

			System.out.println("getMaxXCm() -> "+rslt[0]+"/"+rslt[1]+"/"+rslt[2]);

		}

		return rslt;

	}



	

	/*

	public double [] getMaxXCm( // get fractional center as a "center of mass" inside circle/square from the integer max // USED in lwir

			double [] data,      // rectangular strip of 1/2 of the correlation are with odd rows shifted by 1/2 pixels

			int       center, //  = transform_size - 1;
@@ -2328,7 +2363,7 @@ public class Correlation2d { 
		}

		return rslt;

	}

	

	*/

src/main/java/com/elphel/imagej/tileprocessor/ImageDtt.java

+59 −14
Original line number Diff line number Diff line
@@ -225,20 +225,27 @@ public class ImageDtt extends ImageDttCPU { 
		boolean need_corr = (clt_mismatch != null) || (fcorr_combo_td !=null) || (fcorr_td !=null) ; // (not the only reason)

		// skipping DISPARITY_VARIATIONS_INDEX - it was not used

		if (disparity_map != null){

			for (int i = 0; i<disparity_map.length;i++) if ((disparity_modes & (1 << i)) != 0){

			for (int i = 0; i<disparity_map.length;i++) {

				if (isSliceBit(i) && ((disparity_modes & (1 << i)) != 0)) { 

					if ((i == OVEREXPOSED) && (saturation_imp == null)) {

						continue;

					}

					disparity_map[i] = new double [tilesY*tilesX];

				if ((i >= IMG_TONE_RGB) || ((i >= IMG_DIFF0_INDEX) && (i < (IMG_DIFF0_INDEX + 4)))) {

					need_macro = true;

				}

				if (i <=DISPARITY_STRENGTH_INDEX) {

					if (isCorrBit (i)) {

						need_corr = true;

					}

				} else if (isDiffIndex(i) && needImgDiffs(disparity_modes)){

					disparity_map[i] = new double [tilesY*tilesX];

					need_macro = true;

				} else if (isToneRGBIndex(i) && needTonesRGB(disparity_modes)){

					disparity_map[i] = new double [tilesY*tilesX];

					need_macro = true;

				}

			}

		}



		

		

		if (clt_mismatch != null){

			for (int i = 0; i<clt_mismatch.length;i++){

				clt_mismatch[i] = new double [tilesY*tilesX]; // will use only "center of mass" centers
@@ -373,11 +380,12 @@ public class ImageDtt extends ImageDttCPU { 
				dust_remove,    // boolean   dust_remove,        // Do not reduce average weight when only one image differs much from the average

				false,                 // boolean   calc_textures,

				true);                   // boolean   calc_extra)

			float [][] extra = gpuQuad.getExtra();

			float [][] extra = gpuQuad.getExtra(); // now 4*numSensors

			int num_cams = gpuQuad.getNumCams();

			for (int ncam = 0; ncam < num_cams; ncam++) {

				int indx = ncam + IMG_DIFF0_INDEX;

				if ((disparity_modes & (1 << indx)) != 0){

//				if ((disparity_modes & (1 << indx)) != 0){

				if (needImgDiffs(disparity_modes)){

					disparity_map[indx] = new double [extra[ncam].length];

					for (int i = 0; i < extra[ncam].length; i++) {

						disparity_map[indx][i] = extra[ncam][i];
@@ -386,6 +394,7 @@ public class ImageDtt extends ImageDttCPU { 
			}

			for (int nc = 0; nc < (extra.length - num_cams); nc++) {

				int sindx = nc + num_cams;

				/*

				int indx = nc + IMG_TONE_RGB;

				if ((disparity_modes & (1 << indx)) != 0){

					disparity_map[indx] = new double [extra[sindx].length];
@@ -393,6 +402,17 @@ public class ImageDtt extends ImageDttCPU { 
						disparity_map[indx][i] = extra[sindx][i];

					}

				}

	            */

				int indx = nc + getImgToneRGB(); // IMG_TONE_RGB;

//				if ((disparity_modes & (1 << indx)) != 0){

				if (needTonesRGB(disparity_modes)){

					disparity_map[indx] = new double [extra[sindx].length];

					for (int i = 0; i < extra[sindx].length; i++) {

						disparity_map[indx][i] = extra[sindx][i];

					}

				}



				

			}			

		}

		// does it need non-overlapping texture tiles
@@ -1019,6 +1039,7 @@ public class ImageDtt extends ImageDttCPU { 
		boolean need_macro = false;

		boolean need_corr = (clt_mismatch != null) || (fcorr_combo_td !=null) || (fcorr_td !=null) ; // (not the only reason)

		// skipping DISPARITY_VARIATIONS_INDEX - it was not used

		/*

		if (disparity_map != null){

			for (int i = 0; i<disparity_map.length;i++) if ((disparity_modes & (1 << i)) != 0){

				if ((i == OVEREXPOSED) && (saturation_imp == null)) {
@@ -1033,6 +1054,28 @@ public class ImageDtt extends ImageDttCPU { 
				}

			}

		}

		*/

		if (disparity_map != null){

			for (int i = 0; i<disparity_map.length;i++) {

				if (isSliceBit(i) && ((disparity_modes & (1 << i)) != 0)) { 

					if ((i == OVEREXPOSED) && (saturation_imp == null)) {

						continue;

					}

					disparity_map[i] = new double [tilesY*tilesX];

					if (isCorrBit (i)) {

						need_corr = true;

					}

				} else if (isDiffIndex(i) && needImgDiffs(disparity_modes)){

					disparity_map[i] = new double [tilesY*tilesX];

					need_macro = true;

				} else if (isToneRGBIndex(i) && needTonesRGB(disparity_modes)){

					disparity_map[i] = new double [tilesY*tilesX];

					need_macro = true;

				}

			}

		}



		

		

		if (clt_mismatch != null){

			for (int i = 0; i<clt_mismatch.length;i++){
@@ -1201,7 +1244,8 @@ public class ImageDtt extends ImageDttCPU { 
			int num_cams = gpuQuad.getNumCams();

			for (int ncam = 0; ncam < num_cams; ncam++) {

				int indx = ncam + IMG_DIFF0_INDEX;

				if ((disparity_modes & (1 << indx)) != 0){

//				if ((disparity_modes & (1 << indx)) != 0){

				if (needImgDiffs(disparity_modes)){

					disparity_map[indx] = new double [extra[ncam].length];

					for (int i = 0; i < extra[ncam].length; i++) {

						disparity_map[indx][i] = extra[ncam][i];
@@ -1210,8 +1254,9 @@ public class ImageDtt extends ImageDttCPU { 
			}

			for (int nc = 0; nc < (extra.length - num_cams); nc++) {

				int sindx = nc + num_cams;

				int indx = nc + IMG_TONE_RGB;

				if ((disparity_modes & (1 << indx)) != 0){

				int indx = nc + getImgToneRGB(); // IMG_TONE_RGB;

//				if ((disparity_modes & (1 << indx)) != 0){

				if (needTonesRGB(disparity_modes)){

					disparity_map[indx] = new double [extra[sindx].length];

					for (int i = 0; i < extra[sindx].length; i++) {

						disparity_map[indx][i] = extra[sindx][i];

src/main/java/com/elphel/imagej/tileprocessor/ImageDttCPU.java

+187 −53

File changed.

Preview size limit exceeded, changes collapsed.

src/main/java/com/elphel/imagej/tileprocessor/MacroCorrelation.java

+10 −27
Original line number Diff line number Diff line
@@ -62,6 +62,7 @@ public class MacroCorrelation { 
					mTilesY,   // int tilesY,

					tileSize,  // int tileSize,

					tp.superTileSize, // int superTileSize,

					tp.getNumSensors(),

					tp.isMonochrome(),

					tp.isLwir(),

					tp.isAux(),
@@ -80,7 +81,7 @@ public class MacroCorrelation { 
			final double                             macro_disparity_step,

			final int                                debugLevel){





        int numSensors = geometryCorrection.getNumSensors();

		double [][][] input_data = CLTMacroSetData( // perform single pass according to prepared tiles operations and disparity

				src_scan);           // final CLTPass3d      src_scan, // results of the normal correlations (now expecting infinity)

		if (debugLevel > 0) {
@@ -144,6 +145,8 @@ public class MacroCorrelation { 
		final int mTilesY = (pTilesY + tileSize - 1) / tileSize;

		final int mTiles = mTilesX * mTilesY;

		final int num_chn = 3;

		final int numSensors = tp.getNumSensors();

		final int toneRGB = ImageDtt.getImgToneRGB(numSensors);

  		double     corr_red =          0.5;  // Red to green correlation weight

  		double     corr_blue =         0.2;  // Blue to green correlation weight

  		double [] col_weights = new double[3];
@@ -152,10 +155,10 @@ public class MacroCorrelation { 
		col_weights[1] = corr_blue * col_weights[2];





		final double [][][] input_data = new double [ImageDtt.QUAD][num_chn][mTiles*tileSize*tileSize];

		final int INDX_R0 = ImageDtt.IMG_TONE_RGB;

		final int INDX_B0 = ImageDtt.IMG_TONE_RGB +     ImageDtt.QUAD;

		final int INDX_G0 = ImageDtt.IMG_TONE_RGB + 2 * ImageDtt.QUAD;

		final double [][][] input_data = new double [numSensors][num_chn][mTiles*tileSize*tileSize];

		final int INDX_R0 = toneRGB;

		final int INDX_B0 = toneRGB +     numSensors;

		final int INDX_G0 = toneRGB + 2 * numSensors;



		for (int sub_cam =0; sub_cam < input_data.length; sub_cam++){

			for (int pty = 0; pty < pTilesY; pty++){
@@ -279,34 +282,14 @@ public class MacroCorrelation { 
		}

		double min_corr_selected = clt_parameters.min_corr;



		double [][] disparity_map = new double [ImageDtt.DISPARITY_TITLES.length][]; //[0] -residual disparity, [1] - orthogonal (just for debugging)

//		double [][] disparity_map = new double [ImageDtt.DISPARITY_TITLES.length][]; //[0] -residual disparity, [1] - orthogonal (just for debugging)

		double [][] disparity_map = new double [ImageDtt.getDisparityTitles(geometryCorrection.getNumSensors()).length][]; //[0] -residual disparity, [1] - orthogonal (just for debugging)



		double [][] shiftXY = {{0.0,0.0},{0.0,0.0},{0.0,0.0},{0.0,0.0}};



		double [][][][] clt_corr_combo = null; //    new double [ImageDtt.TCORR_TITLES.length][mTilesY][mTilesX][]; // needed always

		double [][][][][]   clt_corr_partial = null; // [tp.tilesY][tp.tilesX][pair][color][(2*transform_size-1)*(2*transform_size-1)]



		/*

		if (show_corr_partial) {

			clt_corr_partial = new double [mTilesY][mTilesX][][][];

			for (int i = 0; i < mTilesY; i++){

				for (int j = 0; j < mTilesX; j++){

					clt_corr_partial[i][j] = null;

				}

			}

		}



		if (show_corr_combo) {

			clt_corr_combo =  new double [ImageDtt.TCORR_TITLES.length][mTilesY][mTilesX][]; // needed always

			for (int i = 0; i < mTilesY; i++){

				for (int j = 0; j < mTilesX; j++){

					for (int k = 0; k<clt_corr_combo.length; k++){

						clt_corr_combo[k][i][j] = null;

					}

				}

			}

		}

		*/

		

		ImageDtt image_dtt = new ImageDtt(

				geometryCorrection.getNumSensors(),

src/main/java/com/elphel/imagej/tileprocessor/OpticalFlow.java

+5 −5
Original line number Diff line number Diff line
@@ -3227,7 +3227,7 @@ public class OpticalFlow { 
					tilesY,

					true,

					"accumulated_disparity_map-"+nrefine,

					ImageDtt.DISPARITY_TITLES

					ImageDtt.getDisparityTitles(ref_scene.getNumSensors()) // ImageDtt.DISPARITY_TITLES

					);

			// update disparities

			final int disparity_index = ImageDtt.DISPARITY_INDEX_CM; // 2
@@ -3426,7 +3426,7 @@ public class OpticalFlow { 
						tilesY,

						true,

						"accumulated_disparity_map-"+nrefine,

						ImageDtt.DISPARITY_TITLES

						ImageDtt.getDisparityTitles(ref_scene.getNumSensors()) // ImageDtt.DISPARITY_TITLES

						);

			}

			// update disparities
@@ -3637,7 +3637,7 @@ public class OpticalFlow { 
						tilesY,

						true,

						"accumulated_disparity_map-"+nrefine,

						ImageDtt.DISPARITY_TITLES

						ImageDtt.getDisparityTitles(ref_scene.getNumSensors()) // ImageDtt.DISPARITY_TITLES

						);

			}

			// update disparities
@@ -3907,7 +3907,7 @@ public class OpticalFlow { 
				ref_scene.isLwir(),

				clt_parameters.getScaleStrength(ref_scene.isAux()),

				ref_scene.getGPU());

		double[][] disparity_map = new double [ImageDtt.DISPARITY_TITLES.length][];

		double[][] disparity_map = new double [image_dtt.getDisparityTitles().length][];



		int disparity_modes = 

				ImageDtt.BITS_ALL_DISPARITIES |
@@ -4818,7 +4818,7 @@ public class OpticalFlow { 
				ref_scene.isLwir(),

				clt_parameters.getScaleStrength(ref_scene.isAux()),

				ref_scene.getGPU());

		double[][] disparity_map = new double [ImageDtt.DISPARITY_TITLES.length][];

		double[][] disparity_map = new double [image_dtt.getDisparityTitles().length][];



		int disparity_modes = 

				ImageDtt.BITS_ALL_DISPARITIES |