Extracting radial pattern

src/main/java/com/elphel/imagej/orthomosaic/ComboMatch.java

@@ -46,8 +46,8 @@ public class ComboMatch {
 	public static boolean openTestPairGps(
 			 CLTParameters    clt_parameters,
 			 GPUTileProcessor gpu_tile_processor, // initialized by the caller
-			 boolean    extract_mines,
-			 int        debugLevel) {
+			 boolean          extract_objects,
+			 int              debugLevel) {
 		GPU_TILE_PROCESSOR = gpu_tile_processor;
 		PairwiseOrthoMatch pairwiseOrthoMatch = null;
 		String [] pair_names = new String[2];
@@ -55,17 +55,31 @@ public class ComboMatch {
 		String [] gpu_spair = null;
 		boolean use_marked_image =    false; // will be set if found
 		String [] all_scenes = null;
-		if (imp_sel != null) {
-			pair_names = new String[2];
-			String [][] all_scenes1=new String[1][];
-			pairwiseOrthoMatch = getPairFromSelectedImage(
-					imp_sel, // ImagePlus imp,
-					pair_names, // String [] names)
-					all_scenes1); // String [][] all_scenes
-			if (pairwiseOrthoMatch != null) {
-				gpu_spair = pair_names;
-				use_marked_image = true;
-				all_scenes = all_scenes1[0];
+		ArrayList<ObjectLocation> object_list = null;
+		int zool_lev_objects = 0;
+		if (extract_objects) {
+			if (imp_sel == null) {
+				System.out.println("Need multi-sliced image with marked objects as points");
+				return false;
+			} 
+			int [] zoom_lev_objectsp = new int [1];
+			object_list = getObjectsFromSelectedImage(
+					imp_sel,
+					zoom_lev_objectsp); // ImagePlus imp)
+			zool_lev_objects = zoom_lev_objectsp[0];
+		} else {
+			if (imp_sel != null) {
+				pair_names = new String[2];
+				String [][] all_scenes1=new String[1][];
+				pairwiseOrthoMatch = getPairFromSelectedImage(
+						imp_sel, // ImagePlus imp,
+						pair_names, // String [] names)
+						all_scenes1); // String [][] all_scenes
+				if (pairwiseOrthoMatch != null) {
+					gpu_spair = pair_names;
+					use_marked_image = true;
+					all_scenes = all_scenes1[0];
+				}
 			}
 		}
 		// find -L /media/elphel/SSD3-4GB/lwir16-proc/berdich3/linked/linked_1697875868-1697879449-b/ -type f -name "*-GCORR-GEO.tiff" | sort  > GCORR-GEO.list
@@ -254,6 +268,236 @@ public class ComboMatch {
 			maps_collection = new OrthoMapsCollection(files_list_path); // should have ".list" extension
 		}
 		String [] names = maps_collection.getNames();
+		if (object_list != null) {
+			int       corr_size =                128;
+			int       extr_size =     200; // 256;
+			int       zoomout =        2; //  1; //2; // 1; 2 for simulated,1 - for extracted 
+			int       patt_choice =     3;
+			double    phaseCoeff =      0.98;
+			double    min_corr  =       0.01; // 0.0025; // real max >0.005 with scale -500000;
+			int       radius_search =  32;
+			int       radius_centroid = 4;
+			double    lim_rad =        80; // 60; // 40; // 50.0;
+			double    trans_width =    60; // 40; // 20.0;
+			double    reversal_rad =   0; // 26;// 22; // cut at first direction reversal after
+			double    frac_outliers =   0.4;
+			int default_kernel =        1; // 0; // 1; // default kernel choice 25-> 50// 1 for simulated, 0 - for extracted 
+			boolean   only_correlate =  true;
+			int       sub_pattern_index = 0;
+			
+//			String pattern_dir= "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/debug/mines/pattern_25m_zoom1/synthetic/";
+			String pattern_dir= "/media/elphel/NVME/lwir16-proc/ortho_videos/mines_extract/evening_50m_sept12/debug_cross/new_synth/";
+			String [] pattern_files={
+					"patterns_r30.0_e8.0_ir10.0_ie8.0_is-0.1_or45.0_oe30.0_os-1.5_h8_w0.6_s-40.0_200x200.tif",
+					"old_02.tiff",
+					"patterns_r30.0_e24.0_ir10.0_ie8.0_is-0.1_or62.0_oe40.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e8.0_ir10.0_ie8.0_is-0.1_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e8.0_ir10.0_ie8.0_is0.0_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e6.0_ir10.0_ie8.0_is0.0_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e6.0_ir10.0_ie6.0_is-1.0_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e6.0_ir10.0_ie6.0_is-0.8_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e4.0_ir10.0_ie4.0_is-1.0_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e5.0_ir12.0_ie5.0_is-1.0_or45.0_oe30.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e5.0_ir12.0_ie5.0_is-1.0_or45.0_oe20.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e5.0_ir12.0_ie5.0_is-1.0_or50.0_oe40.0_os-1.0_h8_w0.6_s-60.0_200x200.tif", // best
+					"patterns_r30.0_e5.0_ir12.0_ie5.0_is-1.0_or55.0_oe50.0_os-1.0_h8_w0.6_s-60.0_200x200.tif",
+					"patterns_r30.0_e5.0_ir12.0_ie5.0_is0.0_or50.0_oe40.0_os-1.0_h8_w0.6_s-60.0_200x200.tif"
+					};
+			
+			
+			GenericJTabbedDialog gdo = new GenericJTabbedDialog("Set image pair",1200,900);
+			gdo.addNumericField("Correlation size",          corr_size,       0,4,"",
+					"Correlation size, power of 2");
+			gdo.addNumericField("Extracted size",            extr_size,  0,4,"",
+					"size of a square output size");
+//			gdo.addStringField ("Pattern directory",          OrthoMap.pattern_dir, 180, 	"Absolute path including trailing \"/\".");
+//			gdo.addChoice      ("Pattern filename:",          OrthoMap.pattern_files, OrthoMap.pattern_files[patt_choice]);
+			gdo.addStringField ("Pattern directory",          pattern_dir, 180, 	"Absolute path including trailing \"/\".");
+			gdo.addChoice      ("Pattern filename:",          pattern_files, OrthoMap.pattern_files[patt_choice]);
+			gdo.addNumericField("Zoom-out factor",            zoomout,  0,4,"x", "Reduce pattern resolution to match image.");
+			gdo.addNumericField("Phase correlation coefficient",  phaseCoeff,  3,7,"","1.0 - pure phase correlation, 0.0 - regular correlation.");
+			gdo.addNumericField("Minimal correlation",        min_corr,  5,7,"","Minimal correlation value to keep.");
+			gdo.addNumericField("Search radius",              radius_search,    0,4,"pix", "Radius to search for the object center.");
+			gdo.addNumericField("Centroid radius",            radius_centroid,  0,4,"pix", "Limit centroid radius.");
+			gdo.addMessage("Clean-up circular pattern:");
+		
+			gdo.addNumericField("Limit output patter radius", lim_rad,        0,4,"pix", "All zero outside this radius.");
+			gdo.addNumericField("Transition width",           trans_width,    0,4,"pix", "cosine transition width.");
+			gdo.addNumericField("Last reversal radius",       reversal_rad,   0,4,"pix", "Cut at first derivative reversal after this.");
+			gdo.addNumericField("Outlier fraction",           frac_outliers,  3,7,"", "Remove outliers when averaging different directions.");
+			gdo.addCheckbox     ("Only correlate",            only_correlate, "do not generate circular kernel.");
+			gdo.addNumericField("subpattern index",           sub_pattern_index,0,4,"", "0 - full pattern, 1-8 - half-pattern.");
+
+			//only_correlate
+			gdo.showDialog();
+			if (gdo.wasCanceled()) return false;
+			corr_size =                (int) gdo.getNextNumber();
+			extr_size =           (int) gdo.getNextNumber();
+//			OrthoMap.pattern_dir=       gdo.getNextString();		
+//			String pattern_file = OrthoMap.pattern_files[gdo.getNextChoiceIndex()];
+			pattern_dir=                gdo.getNextString();		
+			String pattern_file =       pattern_files[gdo.getNextChoiceIndex()];
+			zoomout=              (int) gdo.getNextNumber();
+			phaseCoeff=                 gdo.getNextNumber();
+			min_corr=                   gdo.getNextNumber();
+			radius_search=        (int) gdo.getNextNumber();
+			radius_centroid=      (int) gdo.getNextNumber();
+			lim_rad=                    gdo.getNextNumber();      
+			trans_width=                gdo.getNextNumber();
+			reversal_rad=               gdo.getNextNumber();
+			frac_outliers=              gdo.getNextNumber();
+			only_correlate =            gdo.getNextBoolean();
+			sub_pattern_index=    (int) gdo.getNextNumber();
+			
+			System.out.println("Will extract objects here, image zoom_level="+zool_lev_objects);
+			for (ObjectLocation ol: object_list) {
+				System.out.println(ol.name+": "+ol.xy_meters[0]+"/"+ol.xy_meters[1]);
+			}
+			maps_collection.reverseRender(
+					object_list, // ArrayList<ObjectLocation> objects,
+					zool_lev_objects);   // int         zoom_level){
+			System.out.println("Object coordinates in source images: ");
+			for (ObjectLocation ol: object_list) {
+				System.out.println(ol.name+": "+ol.getPixels()[0]+"/"+ol.getPixels()[1]);
+			}
+			double [][] object_stack = new double [object_list.size()][];
+			String [] object_titles = new String [object_stack.length];
+			for (int i = 0; i < object_stack.length; i++) {
+				ObjectLocation ol = object_list.get(i);
+				object_titles[i] = ol.getName();
+				object_stack[i] = ol.extractObjectImage(
+						maps_collection,
+						extr_size);
+
+			}
+			ImagePlus imp_obj = ShowDoubleFloatArrays.makeArrays(
+					object_stack,
+					extr_size,
+					extr_size,
+					OrthoMap.removeKnownExtension(imp_sel.getTitle())+"-OBJECT_"+extr_size+"x"+extr_size+".tiff",				
+					object_titles); // test_titles,
+			imp_obj.show();
+			
+			// get pattern(s)
+//			String pattern_path=OrthoMap.pattern_dir+pattern_file;
+			String pattern_path=pattern_dir+pattern_file;
+			ImagePlus imp_pattern = new ImagePlus(pattern_path);
+			int pattern_size =  imp_pattern.getWidth();
+			if (pattern_size == 0) {
+				System.out.println("testPatternCorrelate(): pattern \""+pattern_path+"\" is not found.");
+				return false;
+			}
+			System.out.println("Using pattern file: "+pattern_path+", subpattern="+sub_pattern_index);
+			ImageStack stack_pattern = imp_pattern.getStack();
+			int nSlices = stack_pattern.getSize();
+			double [][] patterns = new double[nSlices][];
+			String [] pattern_labels = new String[nSlices];
+			for (int n = 0; n < patterns.length; n++) {
+				pattern_labels[n]=stack_pattern.getShortSliceLabel(n+1);
+				float [] fpixels_pattern = (float[]) stack_pattern.getPixels(n+1);
+				patterns[n]=new double[fpixels_pattern.length];
+				for (int i = 0; i < fpixels_pattern.length; i++) {
+					patterns[n][i] = fpixels_pattern[i];
+				}
+			}
+			double [] kernel =  OrthoMap.getConvolutionKernel(default_kernel);
+			double [][] centers = new double [object_stack.length][];
+			for (int i = 0; i < object_stack.length; i++) {
+				String    dbg_prefix = "corr_patt_"+i;
+				ObjectLocation ol = object_list.get(i);
+				centers[i] = ObjectLocation.getPatternCenter(
+						object_stack[i], // double [] data,
+						patterns[sub_pattern_index],     // double [] pattern,
+						kernel,          // double [] kernel,
+						zoomout,         // int       zoomout,
+						phaseCoeff,      // double    phaseCoeff,
+						min_corr,        // double    min_corr,
+						radius_search,   // int       radius_search,
+						radius_centroid, // int       radius_centroid,
+						dbg_prefix,      // String    dbg_prefix,
+						debugLevel);     // int       debugLevel)
+				System.out.println(i+": center at "+centers[i][0]+"/"+centers[i][1]+", strength="+centers[i][2]);
+			}
+			PointRoi roi = new PointRoi();
+			roi.setOptions("label");
+			for (int i = 0; i < centers.length; i++) {
+				roi.addPoint(centers[i][0]+extr_size/2, centers[i][1]+extr_size/2, i+1); // ,1);
+			}
+			imp_obj.setRoi(roi);
+			if (only_correlate) {
+				System.out.println("testPatternCorrelate(): correlation DONE, only_correlate is set to true, exiting");
+				return true;
+			}
+			
+//			imp_obj.show();
+			// First process each captured object separately, then compare and average?
+			// TODO: use finer grid, deconvolve with altitude kernel, calculate effective width of the pattern
+			double [][] output_patterns = new double [object_stack.length+1][];
+			String [] output_patt_titles = new String[output_patterns.length];
+			output_patt_titles[output_patt_titles.length-1]= "average";
+			output_patterns[object_stack.length] = new double [corr_size*corr_size];
+//object_titles
+			boolean debug_good = debugLevel > -4;
+			boolean [][] good_pix = debug_good ? (new boolean[object_stack.length][extr_size*extr_size]) : null; // debugging
+			String settings_str = "_fo"+frac_outliers+"_lr"+lim_rad+"_tw"+trans_width+"_rr"+reversal_rad;
+			for (int i = 0; i < object_stack.length; i++) {
+				output_patterns[i] = ObjectLocation.getRadialPattern(
+						object_stack[i], // double [] data,
+						centers[i],      // double [] xy_offs,
+						corr_size,       // int       corr_size, 
+						lim_rad,         // double    lim_rad,     // outside all 0
+						trans_width,     //double    trans_width, // lim_rad-trans_width - start reducing
+						frac_outliers,   // double    frac_outliers,
+						reversal_rad,   // double     reversal_rad,
+						((good_pix != null)? good_pix[i] : null), // boolean [] good,       // null or same size as data
+						debugLevel);     // int       debugLevel)
+				output_patt_titles[i] = object_titles[i];
+				for (int j = 0; j < output_patterns[i].length; j++) {
+					output_patterns[object_stack.length][j] += output_patterns[i][j];
+				}
+			}
+			for (int j = 0; j < output_patterns[0].length; j++) {
+				output_patterns[object_stack.length][j] /= object_stack.length;
+			}
+			ImagePlus imp_out_patt = ShowDoubleFloatArrays.makeArrays(
+					output_patterns,
+					corr_size,
+					corr_size,
+					OrthoMap.removeKnownExtension(imp_sel.getTitle())+"-PATTERN"+settings_str+"_"+corr_size+"x"+corr_size+".tiff",				
+					output_patt_titles); // test_titles,
+			imp_out_patt.show();
+			if (good_pix != null) {
+				double [][] object_stack_masked = new double [2 * object_stack.length][]; 
+				String [] object_titles_masked = new String [2 * object_stack.length];
+				PointRoi roi_masked = new PointRoi();
+				roi_masked.setOptions("label");
+				for (int i = 0; i < centers.length; i++) {
+					roi_masked.addPoint(centers[i][0]+extr_size/2, centers[i][1]+extr_size/2, 2*i+1); // ,1);
+					roi_masked.addPoint(centers[i][0]+extr_size/2, centers[i][1]+extr_size/2, 2*i+2); // ,1);
+				}
+
+				for (int i = 0; i < object_stack.length; i++) {
+					object_stack_masked[2*i + 0] = object_stack[i];
+					object_stack_masked[2*i + 1] = object_stack[i].clone();
+					for (int j = 0; j < object_stack[i].length; j++) {
+						if (!good_pix[i][j]) {
+							object_stack_masked[2*i + 1][j] = Double.NaN;
+						}
+					}
+				}
+				ImagePlus imp_obj_masked = ShowDoubleFloatArrays.makeArrays(
+						object_stack_masked,
+						extr_size,
+						extr_size,
+						OrthoMap.removeKnownExtension(imp_sel.getTitle())+"-OBJECT_MASKED"+settings_str+"_"+extr_size+"x"+extr_size+".tiff",				
+						object_titles_masked); // test_titles,
+				imp_obj_masked.setRoi(roi_masked);
+				imp_obj_masked.show();
+			}
+			return true;
+		}
+		
+		
 		if (use_marked_image) {
 			// verify matching names
 			//all_scenes
@@ -690,6 +934,50 @@ adjusted affines[1] for a pair: 1694564291_293695/1694564778_589341
 		}
 		return null;
 	}
+
+	public static ArrayList<ObjectLocation> getObjectsFromSelectedImage(
+			ImagePlus imp,
+			int [] zoom_levp) {
+		// get zoom from the name and list of the slice names prefixes (before second "_"
+		String title = imp.getTitle();
+		String[] tokens = title.replace(".","_").split("_");
+		for (int i = 0; i < tokens.length; i++) {
+			if (tokens[i].startsWith("zoom")) {
+				int zoom_level = Integer.parseInt(tokens[i].substring(4));
+				System.out.println("zoom level = "+zoom_level);
+				double pix_size = OrthoMap.getPixelSizeMeters (zoom_level); // meters
+				ImageStack stack_scenes = imp.getStack();
+				int nSlices = stack_scenes.getSize();
+				String [] scene_names = new String [nSlices];
+				for (int n = 0; n < scene_names.length; n++) {
+					scene_names[n] = stack_scenes.getSliceLabel(n+1).substring(0,17);
+				}
+				PointRoi pRoi = (PointRoi) imp.getRoi();
+				FloatPolygon fp = pRoi.getContainedFloatPoints();
+				// Use only 2 last markers
+				if (fp == null) {
+					System.out.println("getPairFromSelectedImage(): no markers found");
+					return null;
+				}
+				ArrayList<ObjectLocation> list = new ArrayList<ObjectLocation>();
+				for (int np = 0; np < fp.npoints; np++) {
+					String name = scene_names[pRoi.getPointPosition(np) -1];
+					double [] xy_meters = {
+							fp.xpoints[np] * pix_size, 
+							fp.ypoints[np] * pix_size};
+					list.add(new ObjectLocation (
+							name,        // String name,
+							xy_meters)); // double [] xy_meters)
+				}
+				if (zoom_levp != null) {
+					zoom_levp[0] = zoom_level;
+				}
+				return list;
+			}
+		}
+		return null;
+	}
+	
 	
 	public static String [] getPairChoices(
 			int [][] pairs,
@@ -703,18 +991,6 @@ adjusted affines[1] for a pair: 1694564291_293695/1694564778_589341
 	
 	
 	
-	/*
-	public static TDCorrTile [] rectilinearCorrelate_TD( // scene0/scene1
-			final CLTParameters          clt_parameters,
-			final float  [][]            fpixels, // to check for empty
-			final int                    img_width,
-			Rectangle                    woi, // if null, use full GPU window
-			final double [][][]          affine,  // [2][2][3] affine coefficients to translate common to 2 images
-			final boolean                batch_mode,
-			final int                    debugLevel) {
-	 */
-	
-	
 	public static double [][][] rectilinearVectorField( // scene0/scene1
 			final CLTParameters          clt_parameters,
 			final float  [][]            fpixels, // to check for empty

src/main/java/com/elphel/imagej/orthomosaic/ObjectLocation.java

 package com.elphel.imagej.orthomosaic;
 
+import java.util.ArrayList;
+import java.util.Arrays;
+import java.util.Collections;
+import java.util.Comparator;
+
+import com.elphel.imagej.common.ShowDoubleFloatArrays;
+import com.elphel.imagej.tileprocessor.Correlation2d;
+import com.elphel.imagej.tileprocessor.TileNeibs;
+
+import ij.ImagePlus;
+import ij.gui.PointRoi;
+
 public class ObjectLocation {
 	String name;
 	double [] xy_meters;
+	double [] xy_pixels;
+	
+	public ObjectLocation (
+			String name,
+			double [] xy_meters) {
+		this.name = name;
+		this.xy_meters = xy_meters;
+	}
+	public String getName() {
+		return name;
+	}
+	public double [] getMetric() {
+		return xy_meters;
+	}
+	public int [] getPixels() {
+		return new int[] {
+				(int) Math.round(xy_pixels[0]),
+				(int) Math.round(xy_pixels[1])};
+	};
+	
+	public void setPixels(double [] xy) {
+		xy_pixels = xy;
+	}
+	
+	public double [] extractObjectImage(
+			OrthoMapsCollection maps_collection,
+			int size) {
+		int nmap =maps_collection.getIndex(getName());
+		OrthoMap [] ortho_maps=maps_collection.getMaps();
+		final int width =     ortho_maps[nmap].getImageData().width;
+		final int height =    ortho_maps[nmap].getImageData().height;
+		final float [] src_img =  ortho_maps[nmap].getImageData().data;
+		double [] dcrop = new double [size*size];
+		Arrays.fill(dcrop, Double.NaN);
+		int hsize = size/2;
+		int [] xy = getPixels();
+		boolean has_NaN=false;		
+		for (int y = 0; y < size; y++) {
+			int src_y = y - hsize + xy[1];
+			if ((src_y >= 0) && (y < height)) {
+				for (int x = 0; x < size; x++) {
+					int src_x = x - hsize + xy[0];
+					if ((src_x >= 0) && (x < width)) {
+						double d = src_img[src_x + src_y * width];
+						dcrop[x + size * y] = d;
+						has_NaN |= Double.isNaN(d);
+					} else {
+						has_NaN=true;
+					}
+				}
+			} else {
+				has_NaN=true;
+			}
+		}
+		if (has_NaN) {
+			TileNeibs tn =  new TileNeibs(size,size);
+			OrthoMap.fillNaNs(
+					dcrop, // double [] data,
+					tn, // TileNeibs tn,
+					3); // int min_neibs)
+		}
+		return dcrop;
+	}
+	
+	
+	public static double [] getPatternCenter(
+			double [] data_in,
+			double [] pattern,
+			double [] kernel,
+			int       zoomout,
+			double    phaseCoeff,
+			double    min_corr,
+			int       radius_search,
+			int       radius_centroid,
+			String    dbg_prefix,
+			int       debugLevel){
+		int refine = 1;
+		int extr_size =    (int) Math.sqrt(data_in.length);
+		int corr_size = 1;
+		for (; corr_size <= extr_size; corr_size*=2);
+		corr_size /=2;
+		double [] data = data_in;
+		if (corr_size < extr_size) {
+			data = new double [corr_size*corr_size];
+			int offs = extr_size/2 - corr_size/2;
+			for (int i = 0; i < corr_size; i++) {
+				System.arraycopy(
+						data_in,
+						offs * (extr_size + 1) + i * extr_size,
+						data,
+						i * corr_size,
+						corr_size);
+			}
+		}
+		
+		int pattern_size = (int) Math.sqrt(pattern.length);
+		/*
+		// Already done while extracting
+		boolean has_NaN=false;
+		for (int i = 0; i < data.length; i++) {
+			if (Double.isNaN(data[i])) {
+				has_NaN=true;
+				break;
+			}
+		}
+		if (has_NaN) {
+			TileNeibs tn =  new TileNeibs(corr_size,corr_size);
+			OrthoMap.fillNaNs(
+					data, // double [] data,
+					tn, // TileNeibs tn,
+					3); // int min_neibs)
+		}
+        */
+		if (radius_search > (corr_size/2 - 1)) {
+			radius_search = corr_size/2 - 1;
+			System.out.println("getPatterCenter(): limiting radius_search to "+radius_search);
+		}
+		if (radius_centroid > radius_search) {
+			radius_centroid = radius_search;
+			System.out.println("getPatterCenter(): limiting radius_centroid to "+radius_centroid);
+		}
+		double [] corr_pattern = OrthoMap.patternZoomCropPad(
+				pattern,      // double [] pattern,
+				pattern_size, // int       pattern_size,
+				corr_size,    // int       size,
+				zoomout,      // int       zoomout,
+				false); // true);        // out_normalize); // boolean   normalize) 
+		if (kernel != null) {
+			corr_pattern = OrthoMap.convolveWithKernel(
+					corr_pattern,   // final double [] data,
+					kernel,           // final double [] kernel,
+					corr_size);       // final int width)
+		}
+		double [] corr_out= OrthoMap.correlateWithPattern(
+				data,             // final double [] data,
+				corr_size,        // final int       width,
+				corr_size,        // final int       psize,      // power of 2, such as 64
+				corr_pattern,     // final double [] pattern,    // [psize*psize]
+				false,            // final boolean   convolve,   // convolve, not correlate
+				phaseCoeff,       // final double    phaseCoeff,
+				debugLevel);      // final int       debugLevel)
+		int data_width = 2 *radius_search + 1;
+		int offs = corr_size/2 - radius_search;
+		
+		double [] data_max = new double [data_width*data_width];
+		for (int i = 0; i < data_width; i++) {
+			System.arraycopy(
+					corr_out,
+					offs * (corr_size + 1) + i * corr_size,
+					data_max,
+					i * data_width,
+					data_width);
+		}
+		// TODO: use min_corr and clusterize! Or use fixed fraction (50%) of the maximum ?
+		double [] xys = Correlation2d.getMaxXYCm( // last, average
+				data_max,        // corrs.length-1], // double [] data,
+				data_width,      // int       data_width,      //  = 2 * transform_size - 1;
+				radius_centroid, // double    radius, // 0 - all same weight, > 0 cosine(PI/2*sqrt(dx^2+dy^2)/rad)
+				refine,          // int       refine, //  re-center window around new maximum. 0 -no refines (single-pass)
+				null,            // boolean [] fpn_mask,
+				false,           // boolean    ignore_border, // only if fpn_mask != null - ignore tile if maximum touches fpn_mask
+				false);          // boolean   debug)
+
+		
+		if (dbg_prefix != null) {
+			String []   dbg_titles = {"corr", "orig", "pattern"};
+			double [][] dbg_img = {corr_out, data, corr_pattern};
+			ImagePlus imp_corrs = ShowDoubleFloatArrays.makeArrays(
+					dbg_img,
+					corr_size,
+					corr_size,
+					dbg_prefix+"_pc"+phaseCoeff+"_rad"+radius_centroid,				
+					dbg_titles); // test_titles,
+			PointRoi roi = new PointRoi();
+			roi.setOptions("label");
+			roi.addPoint(xys[0]+corr_size/2, xys[1]+corr_size/2); // ,1);
+			imp_corrs.setRoi(roi);
+			imp_corrs.show();
+		}
+		return xys;
+	}
+
+	public static double [] getRadialPattern(
+			final double [] data,
+			double []  xy_offs,
+			int        corr_size, 
+			double     lim_rad,     // outside all 0
+			double     trans_width, // lim_rad-trans_width - start reducing
+			double     frac_outliers,
+			double     reversal_rad,
+			boolean [] good,       // null or same size as data
+			int        debugLevel) {
+		if (good != null) {
+			Arrays.fill(good, false);
+		}
+		double min_abs_weight = 10.0; // Do not touch very center?
+		int extr_size = (int) Math.sqrt(data.length);
+		double [] rad_patt = new double [corr_size*corr_size];
+		double xc = extr_size/2+xy_offs[0]; 
+		double yc = extr_size/2+xy_offs[1]; 
+		int size_1d = (int) Math.ceil(lim_rad);
+		double lr2 = size_1d * size_1d; 
+		double [] sw =  new double [size_1d + 2];
+		double [] swd = new double [sw.length];
+		double [] swd2 = new double [sw.length];
+//		boolean [] good = new boolean [data.length];
+		double [] distance = new double [data.length];
+		ArrayList<ArrayList<Integer>> rad_lists = new ArrayList<ArrayList<Integer>>();
+		for (int i = 0; i < sw.length; i++) {
+			rad_lists.add(new ArrayList<Integer>());
+		}
+		for (int y = 0; y < extr_size; y++) {
+			double dy = y - yc;
+			double y2 = dy*dy;
+			if (y2 <= lr2) {
+				for (int x = 0; x < extr_size; x++) {
+					double dx = x - xc;
+					double r2 = y2 + dx*dx;
+					if (r2 <= lr2) {
+						int indx = x + extr_size * y;
+						double r = Math.sqrt(r2);
+						int ir = (int) Math.floor(r);
+						distance[indx] = r;
+						rad_lists.get(ir).add(indx);
+						rad_lists.get(ir+1).add(indx);
+						double w1 = r-ir;
+						double w0 = 1.0 - w1;
+						double d = data[indx];
+						sw[ir] +=   w0;
+						swd[ir] +=  w0 * d;
+						swd2[ir] += w0 * d * d;
+						sw[ir+1] +=  w1;
+						swd[ir+1] += w1 * d;
+						swd2[ir+1] += w1 * d * d;
+						if (good != null) {
+							good[indx] = true;
+						}
+					}
+				}
+			}
+		}
+		//size_1d
+		for (int i = 0; i < sw.length; i++) {
+			if (sw[i] > 0) {
+				swd[i] /=  sw[i];
+				swd2[i] /= sw[i];
+			}
+		}
+		double [] sw0 = sw.clone(); // to watch how many were removed (frac_outliers)
+		// now sort each list and later remove from either
+		for (ArrayList<Integer> rad_list:rad_lists) {
+			Collections.sort(rad_list, new Comparator<Integer>() {
+				@Override
+				public int compare(Integer lhs, Integer rhs) {
+					double rhsd = data[rhs];
+					double lhsd = data[lhs];
+					return (rhsd > lhsd) ? -1 : (rhsd < lhsd) ? 1 : 0;
+				}
+			});
+		}
+		
+		// min_weight = 10.0
+		// Starting from very center remove fartherst (first or last) pixel,
+		// until remaining weight drops below min_weight or sw0[]*(1-frac_outliers).
+		// weight of the previous radius may fall also, but likely not too much
+		for (int nr = 1; nr < sw.length; nr++) {
+			double min_weight = Math.max(min_abs_weight, sw0[nr] * (1-frac_outliers));
+			ArrayList<Integer> rad_list = rad_lists.get(nr);
+			ArrayList<Integer> rad_list_prev = (nr > 0)? rad_lists.get(nr-1): null;
+			ArrayList<Integer> rad_list_next = (nr < (sw.length -1))? rad_lists.get(nr+1): null;
+			while (sw[nr] > min_weight) { // will break trying to remove last
+				if (rad_list.size()==0) {
+					System.out.println("BUG! - empty list for nr="+nr);
+					break;
+				}
+				int indx_first = rad_list.get(0);
+				int indx_last = rad_list.get(rad_list.size()-1);
+				boolean remove_first = Math.abs(data[indx_first] - swd[nr]) > Math.abs(data[indx_last] - swd[nr]);
+				int indx_worst = remove_first? indx_first : indx_last;
+				double d = data[indx_worst];
+				double dist = distance[indx_worst];
+				boolean other_next = dist >= nr;
+				double removed_w = other_next ? ( 1 - (dist - nr)) : (1 + (dist - nr));
+				double new_w = sw[nr] - removed_w;
+				if (new_w < min_weight) {
+					break;
+				}
+				swd[nr] =  (swd[nr] *  sw[nr] - removed_w *  d) / new_w;
+				swd2[nr] = (swd2[nr] * sw[nr] - removed_w *  d * d) / new_w;
+				sw[nr] = new_w;
+				rad_list.remove(remove_first ? 0 : (rad_list.size()-1));
+				if (good != null) {
+					good[indx_worst] = false; // seems not used - just for debug/display
+				}
+				// now remove from previous/next list
+				if (dist == nr) {
+					break;
+				}
+				int nr_other = other_next? (nr+1): (nr-1);
+				ArrayList<Integer> rad_list_other = other_next? rad_list_next : rad_list_prev;
+				// check for consistency:
+				int indx_in_other = rad_list_other.indexOf(indx_worst);
+				if (indx_in_other < 0) {
+					System.out.println("BUG in getRadialPattern(), pixel "+indx_worst+" from nr="+nr+" is missing in "+nr_other);
+					break;
+				}
+				double removed_w_other = 1.0 - removed_w;
+				double new_w_other = sw[nr_other] - removed_w_other; // no check for remaining weight here
+				swd[nr_other] =  (swd[nr_other] *  sw[nr_other] - removed_w_other *  d) / new_w_other;
+				swd2[nr_other] = (swd2[nr_other] * sw[nr_other] - removed_w_other *  d * d) / new_w_other;
+				sw[nr_other]= new_w_other;
+				rad_list_other.remove(indx_in_other);
+			}
+		}
+		
+		int lr = size_1d;
+		if (reversal_rad > 0) {
+			int irr = (int) reversal_rad;
+			boolean is_growing = swd[irr+1] > swd[irr];
+			for (lr = irr+1; lr <= size_1d; lr++) {
+				if (is_growing) {
+					if (swd[lr+1] < swd[lr]) {
+						break;
+					}
+				} else {
+					if (swd[lr+1] > swd[lr]) {
+						break;
+					}
+				}
+			}
+		}
+		
+		
+		for (int i = 0; i < lr; i++) {
+			swd[i] -= swd[lr];
+		}
+		for (int i = lr; i <swd.length; i++) {
+			swd[i] = 0;
+		}
+		
+		// TODO: add outlier removal here
+		// apply window
+		if (trans_width > 0) {
+			for (int i = (int) (lim_rad-trans_width); i < sw.length; i++) {
+				if (i < lim_rad) {
+					swd[i] *= 0.5*(1.0 + Math.cos(Math.PI*(i - (lim_rad-trans_width))/trans_width));
+				} else {
+					swd[i] = 0;
+				}
+			}
+		}
+		// generate output pattern
+		int oc = corr_size/2;
+		for (int y = 0; y < corr_size; y++) {
+			double dy = y - oc;
+			double y2 = dy*dy;
+			if (y2 <= lr2) {
+				for (int x = 0; x < corr_size; x++) {
+					double dx = x - oc;
+					double r2 = y2 + dx*dx;
+					if (r2 <= lr2) {
+						double r = Math.sqrt(r2);
+						int ir = (int) Math.floor(r);
+						double w1 = r-ir;
+						double w0 = 1.0 - w1;
+						rad_patt[x + corr_size * y] = w0 * swd[ir] + w1 * swd[ir + 1];
+					}
+				}
+			}
+		}
+		return rad_patt;
+	}
+	
+	
+	
+	
 }

src/main/java/com/elphel/imagej/orthomosaic/OrthoMap.java

@@ -67,6 +67,20 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 			"/media/elphel/NVME/lwir16-proc/ortho_videos/kernel_50_100.tiff"
 	};
 
+	public static String pattern_dir= "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/debug/mines/pattern_25m_zoom1/synthetic/";
+	public static String [] pattern_files=
+		{"patterns_50m_zoom1_200x200.tiff",
+		"patterns_50m_evening_zoom1_200x200.tiff",
+		"patterns_50m_evening_01_zoom1_200x200_00.tiff",
+		"patterns_50m_evening_02_zoom1_200x200_00.tiff",
+		"patterns_50m_evening_03_zoom1_200x200_00.tiff",
+		"patterns_50m_evening_04_zoom1_200x200_00.tiff",
+		"patterns_50m_evening_05_zoom1_200x200_00.tiff",
+		"mine1_zoom0.tiff",
+		"mine2_zoom0.tiff",
+		"mine3_zoom0.tiff",
+		"mine4_zoom0.tiff",
+		"mine6_zoom0.tiff"};
 	
 	
 	public static final String ALT_SUFFIX = "-ALT";
@@ -1095,9 +1109,12 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 	 * @return square kernel as a 1-d double array or null
 	 */
 	public static double [] getConvolutionKernel() {
+		return getConvolutionKernel(kernel_paths.length-1);
+	}
+	public static double [] getConvolutionKernel(int choice) {
 		String kernel_path = null;
 		GenericJTabbedDialog gds = new GenericJTabbedDialog("Select kernel path",1200,400);
-		gds.addChoice("Kernel path:", kernel_paths, kernel_paths[kernel_paths.length-1]);
+		gds.addChoice("Kernel path:", kernel_paths, kernel_paths[choice]);
 		gds.showDialog();
 		if (gds.wasCanceled()) return null;
 		int kernel_index = gds.getNextChoiceIndex();
@@ -3038,6 +3055,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 			double [] wnd_in,
 			double    phaseCoeff) {
 //		phaseCoeff = 0.5;
+		boolean dbg=false;
 		double [] data = data_in.clone();
 		double [] pattern = pattern_in.clone();
 		double [] wnd = wnd_in.clone();
@@ -3047,7 +3065,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 		for (int i = 0; i < data.length; i++) {
 			data[i] *= wnd[i];
 		}
-		{
+		if (dbg) {
 			String [] rslt_titles= {"original","window","windowed","pattern"};
 			ShowDoubleFloatArrays.showArrays(
 					new double[][] {data_orig, wnd, data, pattern},
@@ -3096,7 +3114,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 			}
 		}
 		
-		{
+		if (dbg) {
 			String [] rslt_titles= {"corr", "original","window","windowed","pattern"};
 			ShowDoubleFloatArrays.showArrays(
 					new double[][] {corr_out, data_orig, wnd, data_orig2, pattern_orig},
@@ -3106,7 +3124,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 					convolve?"output_conv":("output_corr"+phaseCoeff),
 					rslt_titles);
 		}
-		System.out.println("testPhaseCorr() done");
+		if (dbg) System.out.println("testPhaseCorr() done");
 		return corr_out;
 	}
 	
@@ -3192,8 +3210,8 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 	
 	public static void testPatternCorrelate(
 			ImagePlus imp_src) {
+		/*
 		String pattern_dir= "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/debug/mines/pattern_25m_zoom1/synthetic/";
-//		String pattern_file= "patterns_50m_zoom1_200x200.tiff";
 		String [] pattern_files=
 			{"patterns_50m_zoom1_200x200.tiff",
 			"patterns_50m_evening_zoom1_200x200.tiff",
@@ -3202,6 +3220,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 			"patterns_50m_evening_03_zoom1_200x200_00.tiff",
 			"patterns_50m_evening_04_zoom1_200x200_00.tiff",
 			"patterns_50m_evening_05_zoom1_200x200_00.tiff"};
+      */			
 		int       zoomout =     2; // 1;
 		int       corr_size = 128; // 256;
 		double    phaseCoeff = 0.5;
@@ -3297,7 +3316,6 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 
 		double [][] corrs_out =     new double[patterns.length][]; 
 		double [][] convolve_out =  new double[patterns.length][]; 
-//		double [][] lim_corr =      new double[patterns.length][];
 		double [][] corr_patterns = new double[patterns.length][];
 		for (int n = 0; n < patterns.length; n++) {
 			corr_patterns[n] = patternZoomCropPad(
@@ -3331,6 +3349,11 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 				height,
 				removeKnownExtension(imp_src.getTitle())+"-PATTERN_CORRS",				
 				patt_titles); // test_titles,
+		
+		
+		
+		
+		
 		if (src_marks != null) {
 			PointRoi roi = new PointRoi();
 			roi.setOptions("label");
@@ -3881,40 +3904,16 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 	public static void testPatternGenerate() {
 		int       half_size = 100;
 		boolean   half_pix = false; // center between pixels
-		boolean   evening_mode = false; // glare around object, use radius 
-		double    radius_out= 42;
-		double    edge_out =  6;
-		double    scale_out = 1.5;
+		boolean   evening_mode = true; // false; // glare around object, use radius 
 		double    radius =   30; // 32; // 32;
 		double    edge=      15; //  4;
 		double    radius_in= 15; // if 0 - skip
 		double    edge_in =   8;
 		double    scale_in = -0.05;
+		double    radius_out= 42;
+		double    edge_out =  6;
+		double    scale_out = 1.5;
 		double    scale =  -200; // black 200
-		if (evening_mode) {
-			radius_out= 42;
-			edge_out =  8;
-			scale_out =-1.5; // relative to scale
-			radius =   30; // 32; // 32;
-			edge=      15; //  4;
-			radius_in= 10; // if 0 - skip
-			edge_in =   8;
-			scale_in = -0.7;
-			scale =   -40; // black 200
-			
-		} else {
-			radius_out= 42;
-			edge_out =  6;
-			scale_out =.08;
-			radius =   30; // 32; // 32;
-			edge=      15; //  4;
-			radius_in= 15; // if 0 - skip
-			edge_in =   8;
-			scale_in = -0.05;
-			scale =  -170; // black 200
-		}
-		
-		
 		
 		boolean   normalize = false;
 		int       halves_number = 8; // hidden mines - number of diameter cuts
@@ -3930,15 +3929,15 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
 		if (gd0.wasCanceled()) return;
 		evening_mode=        gd0.getNextBoolean();
 		if (evening_mode) {
-			radius_out= 62; // 42;
-			edge_out =  40; // 8;
-			scale_out =-1.5; // relative to scale
 			radius =   30; // 32; // 32;
-			edge=       24; // 8; //  4;
+			edge=       8; // 8; //  4;
 			radius_in= 10; // if 0 - skip
 			edge_in =   8;
-			scale_in = -0.1; // 7;
-			scale =   -40; // black 200
+			scale_in = 0; // -0.1; // 7;
+			radius_out= 45; // 42;
+			edge_out =  30; // 8;
+			scale_out =-1.0; // relative to scale
+			scale =   -60; // black 200
 		} else {
 			radius_out= 42;
 			edge_out =  6;

src/main/java/com/elphel/imagej/orthomosaic/OrthoMapsCollection.java

@@ -396,6 +396,64 @@ public class OrthoMapsCollection  implements Serializable{
 		return bounds;
 	}
 	
+	/**
+	 * Use object coordinates (now metric) on the composite image to
+	 * determine corresponding pixel coordinates of the source images
+	 * @param objects ArrayList of objects with the scene name and xy coordinates
+	 *                in the composite image. Will fill in pixel coordinates of the objects 
+	 * @param zoom_level zoom level of the composite image
+	 */
+	public void reverseRender(
+			ArrayList<ObjectLocation> objects,
+			int         zoom_level){
+		int [][] bounds = getBoundsPixels( // should be for rectified, {-bounds[0][0], -bounds[0][1]} - exact center
+				zoom_level,
+				null); 
+		int width =  bounds[0][1] - bounds[0][0]; // bounds[x][0] - negative 
+		int height = bounds[1][1] - bounds[1][0];
+		/*
+		int [] indices = new int [ortho_maps.length]; //maybe will get rid of
+		for (int i = 0; i < indices.length; i++) {
+			indices[i] = i;
+		}
+		*/
+		final int reference_index = 0; 
+		for (ObjectLocation ol: objects) {
+			int nmap =getIndex(ol.getName());
+			final double scale = 1.0/OrthoMap.getPixelSizeMeters(zoom_level);
+			final double src_scale = 1.0/OrthoMap.getPixelSizeMeters(ortho_maps[nmap].orig_zoom_level); // pix per meter
+			double [][] mbounds = ortho_maps[nmap].getBoundsMeters(true); // keep original bounds
+			double [] enu_offset = ortho_maps[reference_index].enuOffsetTo(ortho_maps[nmap]);
+			double [] scaled_out_center = { // xy center to apply affine to
+					-bounds[0][0] + scale * enu_offset[0],
+					-bounds[1][0] - scale * enu_offset[1]};
+			double [] metric_center = { // back to meters to match objects
+					scaled_out_center[0]/scale,
+					scaled_out_center[1]/scale};
+			double dX = ol.getMetric()[0]-metric_center[0]; // matching dX, dY in renderMulti()
+			double dY = ol.getMetric()[1]-metric_center[1];
+			
+//			final int [][] obounds = new int [2][2]; // output (rectified, combined) image bounds, relative to thje top-left
+//			for (int n = 0; n< 2; n++) {
+//				obounds[n][0] = (int) Math.floor(scaled_out_center[n] + scale*mbounds[n][0]);
+//				obounds[n][1] = (int) Math.ceil (scaled_out_center[n] + scale*mbounds[n][1]);
+//			}
+			// Output window size (keep original affine - OK if will not exactly fit)
+//			final int ownd_width =  obounds[0][1] - obounds[0][0];
+//			final int ownd_height = obounds[1][1] - obounds[1][0];
+//			final int ownd_len = ownd_width * ownd_height;
+			double [][] src_bounds=ortho_maps[nmap].getBoundsMeters (true); // using original affines
+			final double [] src_center = {-src_bounds[0][0],-src_bounds[1][0]}; // x,y center offset in the source image
+			final double [][] affine = ortho_maps[nmap].affine; // only here use provided
+//			final int src_width =      ortho_maps[nmap].getImageData().width;
+//			final int src_height =     ortho_maps[nmap].getImageData().height;
+//			final float [] src_img =   ortho_maps[nmap].getImageData().data; // FIXME: will not use
+			double [] xy_src = { // pixels of the source image
+					src_scale * (affine[0][0]*dX + affine[0][1]*dY + affine[0][2] + src_center[0]),
+					src_scale * (affine[1][0]*dX + affine[1][1]*dY + affine[1][2] + src_center[1])};
+			ol.setPixels(xy_src);
+		}
+	}
 	
 	
 	public ImagePlus renderMulti (