Updated/tested with center reference

src/main/java/com/elphel/imagej/cameras/CLTParameters.java

@@ -484,7 +484,7 @@ public class CLTParameters {
 
 	
 	
-
+	public boolean    generate_bg =      true;  // Generate background image for 3dmodels
 	public boolean    show_textures    = true;  // show generated textures
 	public boolean    debug_filters    = false;// show intermediate results of filtering
 	// not used anywhere so far
@@ -1624,6 +1624,7 @@ public class CLTParameters {
 		properties.setProperty(prefix+"gmap_crop_extra",            this.gmap_crop_extra+"");            // int     
 		properties.setProperty(prefix+"gmap_tex_pals",              CLTParameters.arr_to_str(this.gmap_tex_pals)); // int[]
 		
+		properties.setProperty(prefix+"generate_bg",                this.generate_bg+"");
 		properties.setProperty(prefix+"show_textures",              this.show_textures+"");
 		properties.setProperty(prefix+"debug_filters",              this.debug_filters+"");
 
@@ -2639,6 +2640,7 @@ public class CLTParameters {
 		if (properties.getProperty(prefix+"gmap_crop_extra")!=null)      this.gmap_crop_extra=Integer.parseInt(properties.getProperty(prefix+       "gmap_crop_extra"));// int     
 		if (properties.getProperty(prefix+"gmap_tex_pals")!=null)        this.gmap_tex_pals=CLTParameters.str_to_iarr(properties.getProperty(prefix+"gmap_tex_pals"));
 		
+		if (properties.getProperty(prefix+"generate_bg")!=null)                   this.generate_bg=Boolean.parseBoolean(properties.getProperty(prefix+"generate_bg"));
 		if (properties.getProperty(prefix+"show_textures")!=null)                 this.show_textures=Boolean.parseBoolean(properties.getProperty(prefix+"show_textures"));
 		if (properties.getProperty(prefix+"debug_filters")!=null)                 this.debug_filters=Boolean.parseBoolean(properties.getProperty(prefix+"debug_filters"));
 
@@ -3900,6 +3902,7 @@ public class CLTParameters {
 		
 		
 		gd.addMessage     ("Earlier 3D generation parameters");
+		gd.addCheckbox    ("Generate background (infinity) image for 3D (not needed for view down)",                 this.generate_bg);
 		gd.addCheckbox    ("Show generated textures",                                                                this.show_textures);
 		gd.addCheckbox    ("show intermediate results of filtering",                                                 this.debug_filters);
 
@@ -5112,6 +5115,7 @@ public class CLTParameters {
 		
 		this.gmap_tex_pals =           CLTParameters.str_to_iarr(gd.getNextString());
 		
+		this.generate_bg=           gd.getNextBoolean();
 		this.show_textures=         gd.getNextBoolean();
 		this.debug_filters=         gd.getNextBoolean();
 		this.min_smth=              gd.getNextNumber();

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

@@ -237,6 +237,10 @@ public class ErsCorrection extends GeometryCorrection {
 	public double [] getCameraATR() {
 		return camera_atr;
 	}
+	
+	public double [][] getErsXYZATR_dt() {
+		return new double[][] {ers_wxyz_center_dt,ers_watr_center_dt};
+	}
 	public double [] getErsXYZ_dt() {
 		return ers_wxyz_center_dt;
 	}
@@ -378,7 +382,7 @@ public class ErsCorrection extends GeometryCorrection {
 		scenes_poses = new HashMap <String, XyzAtr>();
 	}
 
-	public void removeScene(String timestamp) {
+	public void removeScene(String timestamp) { // not used
 		scenes_poses.remove(timestamp);
 	}
 
@@ -393,20 +397,45 @@ public class ErsCorrection extends GeometryCorrection {
 		return scenes;
 	}
 
-	public void addScene(String timestamp, XyzAtr scene) {
-		scenes_poses.put(timestamp, scene);
+	public void addScene(String timestamp) { // will it work?
+//		scenes_poses.put(timestamp, null);
+		scenes_poses.remove(timestamp);
+	}
+	
+	public void addScene(String timestamp, XyzAtr scene) { // will it work?
+		if (scene != null) {
+			scenes_poses.put(timestamp, scene);
+		} else {
+			scenes_poses.remove(timestamp);
+		}
 	}
+
 	public void addScene(String timestamp, double [] xyz, double [] atr) {
 		scenes_poses.put(timestamp, new XyzAtr(xyz, atr));
 	}
+	public void addScene(String timestamp, double [][] xyzatr) {
+		scenes_poses.put(timestamp, new XyzAtr(xyzatr[0], xyzatr[1]));
+	}
 	
 	public void addScene(String timestamp, double [] xyz, double [] atr, double [] ers_xyz_dt, double [] ers_atr_dt) {
 		scenes_poses.put(timestamp, new XyzAtr(xyz, atr, ers_xyz_dt, ers_atr_dt));
 	}
+	
+	public void addScene(String timestamp, double [][] xyzatr, double [][] ers_xyzatr_dt) {
+		scenes_poses.put(timestamp, new XyzAtr(xyzatr[0], xyzatr[1], ers_xyzatr_dt[0], ers_xyzatr_dt[1]));
+	}
+	
 	//not used
 	public void addScene(String timestamp, double [] xyz, double [] atr, double [] ers_xyz_dt, double [] ers_atr_dt, double [] ers_xyz_d2t, double [] ers_atr_d2t) {
 		scenes_poses.put(timestamp, new XyzAtr(xyz, atr, ers_xyz_dt, ers_atr_dt, ers_xyz_d2t, ers_atr_d2t));
 	}
+
+	public void addScene(String timestamp, double [][] xyzatr, double [][] ers_xyzatr_dt, double [][] ers_xyzatr_d2t) {
+		scenes_poses.put(timestamp, new XyzAtr(
+				xyzatr[0], xyzatr[1],
+				ers_xyzatr_dt[0], ers_xyzatr_dt[1],
+				ers_xyzatr_d2t[0], ers_xyzatr_d2t[1]));
+	}
 	
 	public String getEarliestScene() {
 		return Collections.min(scenes_poses.keySet());
@@ -419,7 +448,10 @@ public class ErsCorrection extends GeometryCorrection {
 	
 	public double[][] getSceneXYZATR(String timestamp) {
 		XyzAtr scene = scenes_poses.get(timestamp);
-		if (scene == null) return null;
+		if (scene == null) {
+//			if timestamp.equals()
+			return null;
+		}
 		return scene.getXYZATR();
 	}
 	public double[] getSceneXYZ(String timestamp) {
@@ -432,6 +464,14 @@ public class ErsCorrection extends GeometryCorrection {
 		if (scene == null) return null;
 		return scene.getATR();
 	}
+
+	public double[][] getSceneErsXYZATR_dt(String timestamp) {
+		XyzAtr scene = scenes_poses.get(timestamp);
+		if (scene == null) return null;
+		return scene.getErsXYZATR_dt();
+	}
+	
+	
 	public double[] getSceneErsXYZ_dt(String timestamp) {
 		XyzAtr scene = scenes_poses.get(timestamp);
 		if (scene == null) return null;

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

@@ -31,6 +31,7 @@ import org.apache.commons.math3.geometry.euclidean.threed.Rotation;
 import com.elphel.imagej.calibration.CalibrationFileManagement;
 import com.elphel.imagej.cameras.CLTParameters;
 import com.elphel.imagej.cameras.ColorProcParameters;
+import com.elphel.imagej.cameras.EyesisCorrectionParameters;
 import com.elphel.imagej.common.DoubleGaussianBlur;
 import com.elphel.imagej.common.ShowDoubleFloatArrays;
 import com.elphel.imagej.gpu.GpuQuad;
@@ -77,10 +78,192 @@ public class Interscene {
 	public Interscene (OpticalFlow opticalFlow) {
 		this.opticalFlow = opticalFlow;
 	}
+/*
+ * Should be done before calling getSceneRange():
+		for (int scene_index =  ref_index - 1; scene_index >= 0 ; scene_index--) {
+			// to include ref scene photometric calibration
+			quadCLTs[scene_index] = quadCLTs[ref_index].spawnNoModelQuadCLT(
+					set_channels[scene_index].set_name,
+					clt_parameters,
+					colorProcParameters, //
+					threadsMax,
+					debugLevel-2);
+		} // split cycles to remove output clutter
 
+ */
+	/**
+	 * Calculate subrange and a center index, starting from the end of the input range
+	 * (former reference). In addition to the IMS data this estimation requires
+	 * knowledge of the ground level. Initial estimate is provided as input,
+	 * then, after knowing the center (reference) index, the DSI will be calculated
+	 * and compared with the estimate. If it differs too much the range will be
+	 * recalculated for the updated ground level. 
+	 * @param clt_parameters  processing parameters
+	 * @param scenes          scenes with IMS data
+	 * @param ground_level    estimated ground level ASL
+	 * @param src_range       a pair of {start, end} scene indices (0<=start,end<
+	 * @param debugLevel      debug level
+	 * @return  {start,end, reference} indices, end index is now the same as src_range[1]          
+	 */
+	public static int []  getSceneRange(
+			final CLTParameters          clt_parameters,
+			final QuadCLT[]              scenes, //
+			final double                 ground_level,
+			final int []                 src_range, // range of scenes to consider {start,end}
+			final int                    debugLevel) {
+		int tilesX =  scenes[src_range[1]].getTileProcessor().getTilesX();
+        int tile_size = scenes[src_range[1]].getTileProcessor().getTileSize();
+		double min_offset =         0.0; // clt_parameters.imp.min_offset;
+		double max_offset =         clt_parameters.imp.max_rel_offset * tilesX * tile_size;
+		double max_roll =           clt_parameters.imp.max_roll_deg*Math.PI/180.0;
+		double max_zoom_diff =      clt_parameters.imp.max_zoom_diff;
+		
+		return null;
+	}
 	
-	public static boolean setInitialOrientationsIms(
+	/**
+	 * Calculate centered initial orientations, return center index
+	 * and also set to config files
+	 * @param clt_parameters
+	 * @param min_num_scenes
+	 * @param colorProcParameters
+	 * @param debayerParameters
+	 * @param rgbParameters
+	 * @param quadCLT_main
+	 * @param quadCLTs
+	 * @param ref_index
+	 * @param set_channels
+	 * @param batch_mode
+	 * @param earliest_scene
+	 * @param start_ref_pointers
+	 * @param threadsMax
+	 * @param updateStatus
+	 * @param debugLevel
+	 * @return
+	 */
+	public static int setInitialOrientationsCenterIms(
+			final CLTParameters                          clt_parameters,
+			int                                          min_num_scenes,
+			final ColorProcParameters                    colorProcParameters,
+			EyesisCorrectionParameters.DebayerParameters debayerParameters,			
+			EyesisCorrectionParameters.RGBParameters     rgbParameters,
+			QuadCLT                                      quadCLT_main, // tiles should be set			
+			final QuadCLT[]                              quadCLTs, //
+			final int                                    ref_index,
+			final QuadCLT.SetChannels []                 set_channels,
+			final boolean                                batch_mode,
+			int                                          earliest_scene,
+			int []                                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
+			final int                                    threadsMax,  // int               threadsMax,
+			final boolean                                updateStatus,
+			final int                                    debugLevel) {
+		// go twice, then inverse 1-st part.
+		// when going first - do not go beyond center (center by IMS coordinates) or just index?
+		// Index is easier, no need to change
+//		QuadCLT[]              quadCLTs_half = new QuadCLT[(quadCLTs.length+1)/2];
+//		QuadCLT.SetChannels [] set_channels_half;
+		int [] start_ref_pointers1 = new int[2];
+		int cent_index =        earliest_scene + (ref_index - earliest_scene) / 2;
+		int min_num_scenes_half =        min_num_scenes / 2; 
+// TODO: add ref scene itself to the list?		
+		cent_index = setInitialOrientationsIms(
+				clt_parameters,      // final CLTParameters          clt_parameters,
+				min_num_scenes_half,     // int                          min_num_scenes,
+				colorProcParameters, // final ColorProcParameters    colorProcParameters,
+				quadCLTs,            // final QuadCLT[]              quadCLTs, //
+				ref_index,           // final int                    ref_index,
+				set_channels,        // final QuadCLT.SetChannels [] set_channels,
+				batch_mode,          // final boolean                batch_mode,
+				cent_index,     // int                          earliest_scene,
+				start_ref_pointers1, // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
+				threadsMax,          // final int                    threadsMax,  // int               threadsMax,
+				updateStatus,        // final boolean                updateStatus,
+				debugLevel);         // final int                    debugLevel)
+		if (cent_index < 0) {
+			System.out.println("setInitialOrientationsIms() first half failed. Consider more graceful bail out.");
+			start_ref_pointers[0] = start_ref_pointers1[0];
+			return cent_index;
+		}
+		// create new dsi for quadCLTs[earliest_scene1]
+		OpticalFlow.buildRefDSI( // returned is a different instance than input -FIXED
+				clt_parameters,                   // CLTParameters                                 clt_parameters,
+				false,                            // boolean  fast,
+				debayerParameters,                // EyesisCorrectionParameters.DebayerParameters  debayerParameters,
+				colorProcParameters,              // ColorProcParameters                           colorProcParameters,
+			    rgbParameters, // EyesisCorrectionParameters.RGBParameters      rgbParameters,
+				batch_mode,                       // boolean                             batch_mode,
+				set_channels[cent_index].set_name,// String                              set_name,
+				quadCLT_main,                     // QuadCLT                             quadCLT_main, // tiles should be set
+				quadCLTs[cent_index],             // QuadCLT                             quadCLT_ref, // tiles should be set
+			    threadsMax,                       // final int                                     threadsMax,  // maximal number of threads to launch
+			    updateStatus,                     // final boolean    updateStatus,
+			    debugLevel);                      // final int        debugLevel);
+// TODO: add {start,end} pointers to quadCLTs[cent_index]
+// TODO: add pointer to center to ref_index
+		
+		
+		
+		int [] start_ref_pointers2 = new int[2];
+		int earliest_scene2 = setInitialOrientationsIms(
+				clt_parameters,      // final CLTParameters          clt_parameters,
+				min_num_scenes_half, // int                          min_num_scenes,
+				colorProcParameters, // final ColorProcParameters    colorProcParameters,
+				quadCLTs,            // final QuadCLT[]              quadCLTs, //
+				cent_index,          // final int                    ref_index,
+				set_channels,        // final QuadCLT.SetChannels [] set_channels,
+				batch_mode,          // final boolean                batch_mode,
+				earliest_scene,      // int                          earliest_scene,
+				start_ref_pointers2, // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
+				threadsMax,          // final int                    threadsMax,  // int               threadsMax,
+				updateStatus,        // final boolean                updateStatus,
+				debugLevel);         // final int                    debugLevel)
+		if (earliest_scene2 < 0) {
+			System.out.println("setInitialOrientationsIms() second half failed. Consider more graceful bail out.");
+			start_ref_pointers[0] = start_ref_pointers1[0];
+			return cent_index;
+		}
+		// invert first half, reference to the cent_index, add to cent_index map, generate ref_index ponter and cent_index,
+		// write config for both ref_index and cent_index scenes
+		ErsCorrection ers_reference = quadCLTs[ref_index].getErsCorrection();
+		ErsCorrection ers_center =    quadCLTs[cent_index].getErsCorrection();
+		String cent_ts = quadCLTs[cent_index].getImageName();
+		double [][] center_xyzatr = new double [][] {ers_reference.getSceneXYZ(cent_ts), ers_reference.getSceneATR(cent_ts)};
+		double [][] inv_cent_xyzatr = ErsCorrection.invertXYZATR(center_xyzatr);
+		for (int scene_index = cent_index; scene_index <= ref_index; scene_index++) { // include cent_index itself to the map
+			double [][] scene_xyzatr,dxyzatr_dt;
+			if (scene_index == ref_index) {
+				scene_xyzatr = new double [2][3];
+				dxyzatr_dt = ers_reference.getErsXYZATR_dt();
+			} else {
+				String ts = quadCLTs[scene_index].getImageName();
+				scene_xyzatr = ers_reference.getSceneXYZATR(ts);
+				dxyzatr_dt =   ers_reference.getSceneErsXYZATR_dt(ts);
+			}
+			double [][] scene_cent_xyzatr = ErsCorrection.combineXYZATR(scene_xyzatr, inv_cent_xyzatr);
+			ers_center.addScene(quadCLTs[scene_index].getImageName(),
+					scene_cent_xyzatr,
+					dxyzatr_dt  // ers_scene.getErsXYZATR_dt(),		
+					);
+		}
+		// set pointers
+		quadCLTs[ref_index].setRefPointer(cent_ts);
+		quadCLTs[cent_index].setFirstLastPointers(quadCLTs[earliest_scene2],quadCLTs[ref_index]);		
+		quadCLTs[ref_index].set_orient(1); // first orientation // applicable to the center?
+		quadCLTs[ref_index].set_accum(0);  // reset accumulations ("build_interscene") number
+		quadCLTs[ref_index].saveInterProperties( // save properties for interscene processing (extrinsics, ers, ...)  // null pointer
+				null, // String path,             // full name with extension or w/o path to use x3d directory
+				debugLevel+1);
+		quadCLTs[cent_index].set_orient(1); // first orientation
+		quadCLTs[cent_index].set_accum(0);  // reset accumulations ("build_interscene") number
+		quadCLTs[cent_index].saveInterProperties( // save properties for interscene processing (extrinsics, ers, ...)  // null pointer
+				null, // String path,             // full name with extension or w/o path to use x3d directory
+				debugLevel+1);
+		return earliest_scene2;
+	}
+	
+	public static int setInitialOrientationsIms(
 			final CLTParameters          clt_parameters,
+			int                          min_num_scenes,
 			final ColorProcParameters    colorProcParameters,
 			final QuadCLT[]              quadCLTs, //
 			final int                    ref_index,
@@ -92,10 +275,10 @@ public class Interscene {
 			final boolean                updateStatus,
 			final int                    debugLevel) {
 		double maximal_series_rms = 0.0;
-		double  min_ref_str =       clt_parameters.imp.min_ref_str;
-		double  min_ref_frac=       clt_parameters.imp.min_ref_frac;		
-    	int min_num_scenes =        clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
-    	int max_num_scenes =        clt_parameters.imp.max_num_scenes; // cut longer series
+		double  min_ref_str =     clt_parameters.imp.min_ref_str;
+		double  min_ref_frac=     clt_parameters.imp.min_ref_frac;		
+//    	int min_num_scenes =      clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
+    	int max_num_scenes =      clt_parameters.imp.max_num_scenes; // cut longer series
     	boolean ims_use =         clt_parameters.imp.ims_use;
 		double [] ims_ortho =     clt_parameters.imp.ims_ortho;
 		double [] ims_mount_atr = clt_parameters.imp.getImsMountATR(); // converts to radians
@@ -106,7 +289,7 @@ public class Interscene {
     	
     	if (!ims_use) {
     		System.out.println("setInitialOrientationsIms(): IMS use is disabled");
-    		return false;
+    		return -1;
     	}
 		if (clt_parameters.ilp.ilma_3d) {
 			System.out.println("*** setInitialOrientationsIms(): clt_parameters.ilp.ilma_3d is TRUE. You may want to disable it when using IMS");
@@ -187,6 +370,11 @@ public class Interscene {
 		double [] ref_abs_atr_enu = Imx5.quatToCamAtr(cam_quat_ref_enu);
 		double [][] ims_ref_xyzatr_enu = {ZERO3, ref_abs_atr_enu};
 		double [][] last_corr_xyzatr = {ZERO3,ZERO3};
+		ers_reference.addScene(quadCLTs[ref_index].getImageName(), // add reference scene (itself) too
+				scenes_xyzatr[ref_index][0],
+				scenes_xyzatr[ref_index][1],
+				ZERO3, // ers_scene.getErsXYZ_dt(),		
+				ZERO3); // ers_scene.getErsATR_dt()		
 		for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
 			if ((ref_index - scene_index) >= max_num_scenes){
 				earliest_scene = scene_index + 1;
@@ -195,7 +383,7 @@ public class Interscene {
 				}
 				// set this and all previous to null
 				for (; scene_index >= 0 ; scene_index--) {
-					ers_reference.addScene(quadCLTs[scene_index].getImageName(), null);
+					ers_reference.addScene(quadCLTs[scene_index].getImageName()); // remove
 				}
 				break;
 			}
@@ -320,7 +508,7 @@ public class Interscene {
 					}
 					// set this and all previous to null
 					for (; scene_index >= 0 ; scene_index--) {
-						ers_reference.addScene(quadCLTs[scene_index].getImageName(), null);
+						ers_reference.addScene(quadCLTs[scene_index].getImageName()); // remove
 						//					quadCLTs[scene_index] = null; // completely remove early scenes?
 					}
 					break;
@@ -377,7 +565,7 @@ public class Interscene {
 			if (start_ref_pointers != null) {
 				start_ref_pointers[0] = earliest_scene;
 			}
-			return false;
+			return -1;
 		}
 		if (earliest_scene > 0) {
 			System.out.println("setInitialOrientationsIms(): Not all scenes matched, earliest useful scene = "+earliest_scene+
@@ -398,7 +586,8 @@ public class Interscene {
     	quadCLTs[ref_index].getErsCorrection().setErsDt( // set for ref also (should be set before non-ref!)
     			dxyzatr_dt[ref_index][0], // (ers_use_xyz? dxyzatr_dt[nscene][0]: ZERO3), //, // dxyzatr_dt[nscene][0], // double []    ers_xyz_dt,
     			dxyzatr_dt[ref_index][1]); // dxyzatr_dt[nscene][1]); // double []    ers_atr_dt)(ers_scene_original_xyz_dt);
-		for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
+//		for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
+		for (int scene_index =  ref_index; scene_index >= earliest_scene ; scene_index--) {
 			ers_reference.addScene(quadCLTs[scene_index].getImageName(),
 					scenes_xyzatr[scene_index][0],
 					scenes_xyzatr[scene_index][1],
@@ -412,13 +601,14 @@ public class Interscene {
 		quadCLTs[ref_index].saveInterProperties( // save properties for interscene processing (extrinsics, ers, ...)  // null pointer
 				null, // String path,             // full name with extension or w/o path to use x3d directory
 				debugLevel+1);
-		return true;
+		return earliest_scene;
 	}
 
 
 	
-	public static boolean setInitialOrientations(
+	public static int setInitialOrientations(
 			final CLTParameters          clt_parameters,
+			int                          min_num_scenes,
 			final ColorProcParameters    colorProcParameters,
 			final QuadCLT[]              quadCLTs, //
 			final int                    ref_index,
@@ -437,7 +627,7 @@ public class Interscene {
 		double maximal_series_rms = 0.0;
 		double  min_ref_str =       clt_parameters.imp.min_ref_str;
 		double  min_ref_frac=       clt_parameters.imp.min_ref_frac;		
-    	int min_num_scenes =        clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
+//    	int min_num_scenes =        clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
     	int max_num_scenes =        clt_parameters.imp.max_num_scenes; // cut longer series
     	
     	double boost_max_short = 2.0; // 
@@ -512,6 +702,11 @@ public class Interscene {
 		boolean after_spiral = false;
 		boolean got_spiral = false;
 		int    search_rad =        clt_parameters.imp.search_rad;    // 10;
+		ers_reference.addScene(quadCLTs[ref_index].getImageName(), // add reference scene (itself) too
+				scenes_xyzatr[ref_index][0],
+				scenes_xyzatr[ref_index][1],
+				ZERO3, // ers_scene.getErsXYZ_dt(),		
+				ZERO3); // ers_scene.getErsATR_dt()		
 		for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
 			if ((ref_index - scene_index) >= max_num_scenes){
 				earliest_scene = scene_index + 1;
@@ -520,7 +715,7 @@ public class Interscene {
 				}
 				// set this and all previous to null
 				for (; scene_index >= 0 ; scene_index--) {
-					ers_reference.addScene(quadCLTs[scene_index].getImageName(), null);
+					ers_reference.addScene(quadCLTs[scene_index].getImageName()); // remove
 				}
 				break;
 			}
@@ -555,9 +750,9 @@ public class Interscene {
 					// set this and all previous to null
 					for (; scene_index >= 0 ; scene_index--) {
 						if (quadCLTs[scene_index] == null) {
-							System.out.println("setInitialOrientations(): quadCLTs["+scene_index+"] is alreadu null!");
+							System.out.println("setInitialOrientations(): quadCLTs["+scene_index+"] is already null!");
 						} else {
-							ers_reference.addScene(quadCLTs[scene_index].getImageName(), null); // null pointer
+							ers_reference.addScene(quadCLTs[scene_index].getImageName()); // remove
 						}
 					}
 					break; // failed with even first before reference
@@ -686,7 +881,7 @@ public class Interscene {
 					}
 					// set this and all previous to null
 					for (; scene_index >= 0 ; scene_index--) {
-						ers_reference.addScene(quadCLTs[scene_index].getImageName(), null);
+						ers_reference.addScene(quadCLTs[scene_index].getImageName()); // remove
 						//					quadCLTs[scene_index] = null; // completely remove early scenes?
 					}
 					break;
@@ -750,7 +945,7 @@ public class Interscene {
 			if (start_ref_pointers != null) {
 				start_ref_pointers[0] = earliest_scene;
 			}
-			return false;
+			return -1;
 		}
 		if (earliest_scene > 0) {
 			System.out.println("setInitialOrientations(): Not all scenes matched, earliest useful scene = "+earliest_scene+
@@ -764,7 +959,7 @@ public class Interscene {
 		quadCLTs[ref_index].saveInterProperties( // save properties for interscene processing (extrinsics, ers, ...)  // null pointer
 				null, // String path,             // full name with extension or w/o path to use x3d directory
 				debugLevel+1);
-		return true;
+		return earliest_scene;
 	}
 	
 	
@@ -773,6 +968,8 @@ public class Interscene {
 			double         mb_max_gain,
 			boolean []     reliable_ref, // null or bitmask of reliable reference tiles			
     		QuadCLT []     quadCLTs,
+    		int            ref_index,
+    		int     []     range,
     		boolean        test_motion_blur,
 			int            debugLevel)
 	{ 
@@ -809,10 +1006,11 @@ public class Interscene {
 		int     margin =      clt_parameters.imp.margin;
         double second_margin = 1.1; // allow slightly larger offset for readjsutments (to be limited by the first pass)
         double boost_max_short = 5.0; // 2.0; // Some bug?
-		int earliest_scene = 0;
+		int earliest_scene = range[0];
+		int last_scene =     range[1];
 	    boolean use_combo_dsi =        clt_parameters.imp.use_combo_dsi;
 	    boolean use_lma_dsi =          clt_parameters.imp.use_lma_dsi;
-    	int ref_index = quadCLTs.length-1;
+	    
 		int tilesX =  quadCLTs[ref_index].getTileProcessor().getTilesX();
         int tilesY =  quadCLTs[ref_index].getTileProcessor().getTilesY();
         int tile_size = quadCLTs[ref_index].getTileProcessor().getTileSize();
@@ -883,14 +1081,15 @@ public class Interscene {
 		
 		scenes_xyzatr[ref_index] = new double[2][3]; // all zeros
 		// should have at least next or previous non-null
-		int debug_scene = -8;
+		int debug_scene = 68; // -8;
 		double maximal_series_rms = 0.0;
 		double [][] mb_vectors_ref = null;
 		TpTask[][] tp_tasks_ref = null;
     	if (debug_ers) {
     		System.out.println("ERS velocities scale mode = '"+dbg_ers_string+"'");
     	}
-        for (int nscene = ref_index; nscene >= earliest_scene; nscene--) {
+    	
+        for (int nscene = last_scene; nscene >= earliest_scene; nscene--) {
 			if (nscene == debug_scene) {
 				System.out.println("nscene = "+nscene);
 				System.out.println("nscene = "+nscene);
@@ -919,13 +1118,13 @@ public class Interscene {
         		quadCLTs,        // QuadCLT []     scenes, // ordered by increasing timestamps
         		ref_index,
         		earliest_scene,  // int            start_scene,
-        		ref_index,       // int            end1_scene,
+        		last_scene,       // int            end1_scene,
         		scenes_xyzatr,   // double [][][]  scenes_xyzatr,
         		half_run_range); // double         half_run_range
         
     	if (debug_ers) {
     		System.out.println();
-    		for (int nscene = ref_index; nscene >= earliest_scene; nscene--) {
+    		for (int nscene = last_scene; nscene >= earliest_scene; nscene--) {
     		System.out.println(String.format("%3d xyz_dt: %9.5f %9.5f %9.5f   atr_dt: %9.5f %9.5f %9.5f",
     				nscene,
     				dxyzatr_dt[nscene][0][0], dxyzatr_dt[nscene][0][1], dxyzatr_dt[nscene][0][2],
@@ -946,7 +1145,20 @@ public class Interscene {
     	
     	boolean [] failed_scenes = new boolean[quadCLTs.length];
     	int num_failed =0;
-        for (int nscene = ref_index; nscene >= earliest_scene; nscene--) {
+    	int [] scene_seq = new int [last_scene-earliest_scene+1];
+    	{
+    		int indx=0;
+    		for (int i = ref_index; i >= earliest_scene; i--) {
+    			scene_seq[indx++] = i;
+    		}
+    		for (int i = ref_index+1; i <= last_scene; i++) {
+    			scene_seq[indx++] = i;
+    		}
+    	}
+    	
+    	
+//        for (int nscene = ref_index; nscene >= earliest_scene; nscene--) {
+       	for (int nscene:scene_seq) {
 			if (nscene == debug_scene) {
 				System.out.println("nscene = "+nscene);
 				System.out.println("nscene = "+nscene);
@@ -1077,7 +1289,15 @@ public class Interscene {
 						mb_max_gain,         // double             mb_max_gain, // 5.0;   // motion blur maximal gain (if more - move second point more than a pixel
 						mb_vectors_ref,      // double [][]        mb_vectors,  // now [2][ntiles];
 						debugLevel);         // int         debugLevel)
-			} else {
+				ers_reference.addScene(ts,
+						scenes_xyzatr[nscene][0], // ref_index
+						scenes_xyzatr[nscene][1],
+						dxyzatr_dt[nscene][0],
+						dxyzatr_dt[nscene][1]
+//						quadCLTs[nscene].getErsCorrection().getErsXYZ_dt(),	// same as dxyzatr_dt[nscene][0], just keep for future adjustments?	
+//						quadCLTs[nscene].getErsCorrection().getErsATR_dt()	// same as dxyzatr_dt[nscene][1], just keep for future adjustments?		
+						);
+			} else { // if (nscene == ref_index)
 				double [][]        mb_vectors = null;
 				scene_xyz_pre = ers_reference.getSceneXYZ(ts);
 				scene_atr_pre = ers_reference.getSceneATR(ts);
@@ -1087,7 +1307,7 @@ public class Interscene {
 				int n_pre =  nscene + avg_len;
 				int n_post = nscene - avg_len;
 				if (n_post < earliest_scene) n_post = earliest_scene;
-				if (n_pre > ref_index) n_pre = ref_index;
+				if (n_pre > last_scene) n_pre = last_scene;
 				double s0=0, sx=0, sx2=0;
 				double [][] sy = new double[2][3], sxy = new double [2][3];
 				for (int n_other = n_post; n_other <= n_pre; n_other++) {
@@ -1165,10 +1385,6 @@ public class Interscene {
 						mb_vectors,  // double [][]    mb_vectors,  // now [2][ntiles];
 						clt_parameters.imp.debug_level);                // 1); // -1); // int debug_level);
 				boolean adjust_OK = scenes_xyzatr[nscene] != null;
-//				if (adjust_OK && (Math.abs(scenes_xyzatr[nscene][1][2]) > max_roll)) {
-//					fail_reason[0] = FAIL_REASON_ROLL;
-//					adjust_OK = false;
-//				}
 				if (adjust_OK && fail_on_zoom_roll) { // check only for initial orientation, do not check on readjustments
 					if (Math.abs(scenes_xyzatr[nscene][1][2]) > max_roll) {
 						fail_reason[0] = FAIL_REASON_ROLL;
@@ -1214,7 +1430,7 @@ public class Interscene {
 						}
 						// set this and all previous to null
 						for (; nscene >= 0 ; nscene--) {
-							ers_reference.addScene(quadCLTs[nscene].getImageName(), null);
+							ers_reference.addScene(quadCLTs[nscene].getImageName()); // remove
 							quadCLTs[nscene] = null; // completely remove early scenes?
 						}
 						break;
@@ -1236,11 +1452,10 @@ public class Interscene {
 							quadCLTs[ref_index].getImageName() + "/" + ts+
 							" Done. RMS="+lma_rms[0]+", maximal so far was "+maximal_series_rms);
 				}
-			}
+			} // if (nscene == ref_index) else
         } // for (int nscene = ref_index; nscene > earliest_scene; nscene--) {
 
         // TODO: after all scenes done, see if any of scenes_xyzatr[] is null, and if fpn_rematch - rematch
-        
 		if (dbg_mb_img != null) {
 			String [] dbg_mb_titles = new String[quadCLTs.length];
 			for (int i = 0; i < quadCLTs.length; i++) if (quadCLTs[i] != null) {
@@ -3335,7 +3550,7 @@ public class Interscene {
 	 * @param selection optional selection to ignore unselected tiles)
 	 * @param margin do not use tiles with centers closer than this to the edges. Measured in pixels.
 	 * @param mb_tau Sensor time constant in seconds (only needed if mb_vectors != null)
-	 * @param mb_max_gain maximal gain fro MB correction (higher MB corrected by increasing offset)
+	 * @param mb_max_gain maximal gain for MB correction (higher MB corrected by increasing offset)
 	 * @param mb_vectors [2][tiles] {dx/dt[tiles],dx/dt[tiles]} motion blur vectors or null 
 	 * @param debug_level debug level
 	 * @return TpTask[][] arrays used to program GPU. With no MB has only one element, with MB - two
@@ -3365,10 +3580,10 @@ public class Interscene {
 		}
 		final double disparity_corr = clt_parameters.imp.disparity_corr; // 04/07/2023 // 0.0; // (z_correction == 0) ? 0.0 : geometryCorrection.getDisparityFromZ(1.0/z_correction);
 //ref_disparity	
-		if (ref_disparity == null) {
-			ref_disparity = ref_scene.getDLS()[use_lma_dsi?1:0];
-		}
 		if (ref_pXpYD == null) {
+			if (ref_disparity == null) {
+				ref_disparity = ref_scene.getDLS()[use_lma_dsi?1:0];
+			}
 			ref_pXpYD = OpticalFlow.transformToScenePxPyD( // full size - [tilesX*tilesY], some nulls
 					null, // final Rectangle [] extra_woi,    // show larger than sensor WOI (or null)
 					ref_disparity, // dls[0],  // final double []   disparity_ref, // invalid tiles - NaN in disparity (maybe it should not be masked by margins?)

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

@@ -763,7 +763,7 @@ public class IntersceneLma {
 		final double [] full_vector = getFullVector(vector);
 		final ErsCorrection ers_ref =   reference_QuadClt.getErsCorrection();
 		final ErsCorrection ers_scene = scene_QuadClt.getErsCorrection();
-		final double [] scene_xyz =     new double[3];;
+		final double [] scene_xyz =     new double[3];
 		final double [] scene_atr  =    new double[3];
 		final double [] reference_xyz = new double[3]; // will stay 0
 		final double [] reference_atr = new double[3]; // will stay 0
@@ -786,6 +786,10 @@ public class IntersceneLma {
 			scene_atr[i] =     full_vector[ErsCorrection.DP_DSAZ +  i];
 			scene_xyz[i] =     full_vector[ErsCorrection.DP_DSX +  i];
 		}
+		if (Double.isNaN(scene_xyz[0])) {
+			System.out.println("getFxDerivs(): scene_xyz[0]=NaN");
+			System.out.println();
+		}
 		ers_scene.setupERS();
 		ers_ref.setupERS();
 		final Matrix [] reference_matrices_inverse = ErsCorrection.getInterRotDeriveMatrices(

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

@@ -52,6 +52,7 @@ public class IntersceneMatchParameters {
 				ims_mount_atr[1] * Math.PI/180,
 				ims_mount_atr[2] * Math.PI/180};
 	}
+	public  boolean center_reference =   false;
 	public  boolean force_ref_dsi =      false; // true;
 	public  boolean force_orientations = false;
 	public  boolean run_ly =             false; // will return just after LY adjustments, skipping all output generation
@@ -161,7 +162,8 @@ public class IntersceneMatchParameters {
 	//------------	
 	
 	// Sky detection/filtering
-	public  boolean sky_recalc =                    false; // force blue sky recalculation even if it exists
+	public  boolean sky_extract =                    true;  // do not try to extract blue sky
+	public  boolean sky_recalc =                     false; // force blue sky recalculation even if it exists
 	public  double  sky_highest_min =                -50; // lowest absolute value should not be higher (requires photometric) 
 	public  double  cold_frac =                        0.005; // this and lower will scale fom by  cold_scale
 	public  double  hot_frac =                         0.9;    // this and above will scale fom by 1.0
@@ -444,6 +446,8 @@ public class IntersceneMatchParameters {
 		gd.addTab         ("Scene Series", "Processing series of scenes and multi-series sets");
 		gd.addMessage  ("Build series options");
 		
+		gd.addCheckbox ("Use center scene as a reference",           this.center_reference,
+				"True for mapping with post-processing. If false, use the last scene as a reference (driving).");
 		gd.addCheckbox ("Force reference scene DSI calculation",     this.force_ref_dsi,
 				"Calculate reference scene DSI even if the file exists.");
 		gd.addCheckbox ("Force egomotion calculation",               this.force_orientations,
@@ -657,6 +661,8 @@ public class IntersceneMatchParameters {
 		
 		gd.addTab("Sky","Featureless sky areas detection and filtering (LWIR only)");
 		gd.addMessage  ("Temperature-related filtering (requires photometric calibration)");
+		gd.addCheckbox ("Try blue sky calculation/use",              this.sky_extract,
+				"Disable if no sky will appear in the image (such as for down-pointed camera).");
 		gd.addCheckbox ("Force blue sky recalculation",              this.sky_recalc,
 				"Recalculate blue sky even if it already exists.");
 		gd.addNumericField("High limit of the scene coldest tile",   this.sky_highest_min, 5,7,"",
@@ -1108,6 +1114,7 @@ public class IntersceneMatchParameters {
 		this.ims_mount_xyz = IntersceneMatchParameters. StringToDoubles(gd.getNextString(), 3);
 		this.sfm_use =                        gd.getNextBoolean();
 		this.sfm_min_base =                   gd.getNextNumber();
+		this.center_reference =               gd.getNextBoolean();
 		this.force_ref_dsi =                  gd.getNextBoolean();
 		this.force_orientations =             gd.getNextBoolean();
 		this.run_ly =                         gd.getNextBoolean();
@@ -1203,6 +1210,7 @@ public class IntersceneMatchParameters {
 		this.strong_strength =                gd.getNextNumber();
 		this.weak_strength =                  gd.getNextNumber();
 		
+		this.sky_extract =              gd.getNextBoolean();
 		this.sky_recalc =               gd.getNextBoolean();
 		this.sky_highest_min =          gd.getNextNumber(); 	
 		this.cold_frac =                gd.getNextNumber();       	
@@ -1461,6 +1469,7 @@ public class IntersceneMatchParameters {
 		properties.setProperty(prefix+"ims_mount_xyz",                 IntersceneMatchParameters.doublesToString(this.ims_mount_xyz));
 		properties.setProperty(prefix+"sfm_use",                       this.sfm_use + "");                     // boolean
 		properties.setProperty(prefix+"sfm_min_base",                  this.sfm_min_base+"");                  // double
+		properties.setProperty(prefix+"center_reference",              this.center_reference + "");               // boolean
 		properties.setProperty(prefix+"force_ref_dsi",                 this.force_ref_dsi + "");               // boolean
 		properties.setProperty(prefix+"force_orientations",            this.force_orientations + "");          // boolean
 		properties.setProperty(prefix+"run_ly",                        this.run_ly + "");                      // boolean
@@ -1563,6 +1572,7 @@ public class IntersceneMatchParameters {
 		properties.setProperty(prefix+"strong_strength",               this.strong_strength+"");               // double
 		properties.setProperty(prefix+"weak_strength",                 this.weak_strength+"");                 // double
 
+		properties.setProperty(prefix+"sky_extract",          this.sky_extract+"");         // boolean
 		properties.setProperty(prefix+"sky_recalc",           this.sky_recalc+"");          // boolean
 		properties.setProperty(prefix+"sky_highest_min",      this.sky_highest_min+"");     // double
 		properties.setProperty(prefix+"cold_frac",            this.cold_frac+"");           // double
@@ -1774,6 +1784,7 @@ public class IntersceneMatchParameters {
 		if (properties.getProperty(prefix+"ims_mount_xyz")!=null) 	   this.ims_mount_xyz= IntersceneMatchParameters.StringToDoubles(properties.getProperty(prefix+"ims_mount_xyz"),3);
 		if (properties.getProperty(prefix+"sfm_use")!=null)            this.sfm_use=Boolean.parseBoolean(properties.getProperty(prefix+"sfm_use"));		
 		if (properties.getProperty(prefix+"sfm_min_base")!=null)       this.sfm_min_base=Double.parseDouble(properties.getProperty(prefix+"sfm_min_base"));
+		if (properties.getProperty(prefix+"center_reference")!=null)              this.center_reference=Boolean.parseBoolean(properties.getProperty(prefix+"center_reference"));		
 		if (properties.getProperty(prefix+"force_ref_dsi")!=null)                 this.force_ref_dsi=Boolean.parseBoolean(properties.getProperty(prefix+"force_ref_dsi"));		
 		if (properties.getProperty(prefix+"force_orientations")!=null)            this.force_orientations=Boolean.parseBoolean(properties.getProperty(prefix+"force_orientations"));		
 		if (properties.getProperty(prefix+"run_ly")!=null)                        this.run_ly=Boolean.parseBoolean(properties.getProperty(prefix+"run_ly"));		
@@ -1881,6 +1892,7 @@ public class IntersceneMatchParameters {
 		if (properties.getProperty(prefix+"strong_strength")!=null)               this.strong_strength=Double.parseDouble(properties.getProperty(prefix+"strong_strength"));
 		if (properties.getProperty(prefix+"weak_strength")!=null)                 this.weak_strength=Double.parseDouble(properties.getProperty(prefix+"weak_strength"));
 		
+		if (properties.getProperty(prefix+"sky_extract")!=null)          this.sky_extract=Boolean.parseBoolean(properties.getProperty(prefix+"sky_extract"));		
 		if (properties.getProperty(prefix+"sky_recalc")!=null)           this.sky_recalc=Boolean.parseBoolean(properties.getProperty(prefix+"sky_recalc"));		
 		if (properties.getProperty(prefix+"sky_highest_min")!=null)      this.sky_highest_min=Double.parseDouble(properties.getProperty(prefix+"sky_highest_min"));
 		if (properties.getProperty(prefix+"cold_frac")!=null)            this.cold_frac=Double.parseDouble(properties.getProperty(prefix+"cold_frac"));
@@ -2114,6 +2126,7 @@ public class IntersceneMatchParameters {
  		imp.ims_mount_xyz                 = this.ims_mount_xyz.clone();
  		imp.sfm_use                       = this.sfm_use;
  		imp.sfm_min_base                  = this.sfm_min_base;
+		imp.center_reference              = this.center_reference;
 		imp.force_ref_dsi                 = this.force_ref_dsi;
 		imp.force_orientations            = this.force_orientations;
 		imp.run_ly                        = this.run_ly;
@@ -2212,6 +2225,7 @@ public class IntersceneMatchParameters {
 		imp.strong_strength               = this.strong_strength;
 		imp.weak_strength                 = this.weak_strength;
 		
+		imp.sky_extract =           this.sky_extract;
 		imp.sky_recalc =            this.sky_recalc;
 		imp.sky_highest_min =       this.sky_highest_min;
 		imp.cold_frac =             this.cold_frac;

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

@@ -3941,6 +3941,7 @@ public class OpticalFlow {
 	
 	public static void buildRefDSI(
 			CLTParameters                                 clt_parameters,
+			boolean                                       fast,
 			EyesisCorrectionParameters.DebayerParameters  debayerParameters,
 			ColorProcParameters                           colorProcParameters,
 			EyesisCorrectionParameters.RGBParameters      rgbParameters,
@@ -3952,8 +3953,8 @@ public class OpticalFlow {
 			final boolean                                 updateStatus,
 			final int                                     debugLevel) {
 		final int        debugLevelInner=clt_parameters.batch_run? -2: debugLevel; // copied from TQ
-    	boolean  photo_each =                clt_parameters.photo_each; //        true;  // perform photogrammetric calibration to equalize pixel values
-
+    	boolean  photo_each =        clt_parameters.photo_each; //        true;  // perform photogrammetric calibration to equalize pixel values
+    	boolean sky_extract =        clt_parameters.imp.sky_extract;
  		double sky_highest_min =     clt_parameters.imp.sky_highest_min;
  		double cold_frac =           clt_parameters.imp.cold_frac;
  		double hot_frac =            clt_parameters.imp.hot_frac;
@@ -4016,9 +4017,30 @@ public class OpticalFlow {
 
 		QuadCLT   dbg_scene = clt_parameters.imp.save_debug_images? quadCLT_ref: null; // use to save debug images if not null
 
-		boolean [] ref_blue_sky = null; // turn off "lma" in the ML output  
+		boolean [] ref_blue_sky = null; // turn off "lma" in the ML output
+		
+		int     gr_max_clust_radius = clt_parameters.gr_max_clust_radius;			  
+		double  disp_scan_start = clt_parameters.disp_scan_start;
+		double  disp_scan_step =  clt_parameters.disp_scan_step;
+		int  disp_scan_count =    clt_parameters.disp_scan_count;
+		boolean no_bg_generate = !clt_parameters.generate_bg;
+		boolean no_lma =         false;			  
+		if (fast) {
+			no_lma =         true;
+			no_bg_generate = true;
+			photo_each =     false;
+			sky_extract =    false;
+		}
+		
 		quadCLT_ref.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
 				clt_parameters,
+				  // adding these parameters for more flexibility in accuracy/speed
+				gr_max_clust_radius, // int      gr_max_clust_radius,			  
+				disp_scan_start,     // double  disp_scan_start,
+				disp_scan_step,      // double  disp_scan_step,
+				disp_scan_count,     // int  disp_scan_count,
+				no_bg_generate,      // boolean no_bg_generate,
+				no_lma,              // boolean no_lma,			  
 				debayerParameters,
 				colorProcParameters,
 				rgbParameters,
@@ -4061,6 +4083,7 @@ public class OpticalFlow {
 			dsi[TwoQuadCLT.DSI_STRENGTH_AUX] =      aux_last_scan[1];
 		}
 		boolean ran_photo_each = false;
+		quadCLT_ref.tp.resetCLTPasses();
 		// perform photometric here, after first DSI
 		if (photo_each && (!quadCLT_ref.isPhotometricThis() || !batch_mode)) {
 //			if (debugLevel > -3) {
@@ -4084,68 +4107,73 @@ public class OpticalFlow {
 			System.out.println("(Re)using photometric calibration from this sequence reference "+quadCLT_ref.getPhotometricScene());
 			quadCLT_ref.setQuadClt(); // just in case ?
 		}
-		quadCLT_ref.setBlueSky  ( // initial BS from single scene
-				max_disparity,
-				max_disparity_strength,
-				sky_seed,           // double sky_seed, //  =       7.0;  // start with product of strength by diff_second below this
-				lma_seed,
-				seed_temp,     //double seed_temp, //         0.5;  // seed colder that this point between min and max temp						
-				sky_lim,            // double sky_lim, //   =      15.0; // then expand to product of strength by diff_second below this
-				lim_temp, //						double lim_temp, //          0.5;  // sky colder that this point between min and max temp						
-				sky_shrink,         // int    sky_shrink, //  =       4;
-				sky_expand_extra,   // int    sky_expand_extra, //  = 100; // 1?
-				sky_bottleneck,     //int    sky_bottleneck, // 
-				sky_reexpand_extra,   //int    sky_reexpand_extra,  // 9; re-expand after bottleneck in addition to how it was shrank
-				cold_scale, // =       0.2;  // <=1.0. 1.0 - disables temperature dependence
-				cold_frac, // =        0.005; // this and lower will scale fom by  cold_scale
-				hot_frac, // =         0.9;    // this and above will scale fom by 1.0
-				min_strength, // =     0.08;
-				seed_rows, // =        5; // sky should appear in this top rows 
-				lowest_sky_row,        //  =   50;// appears that low - invalid, remove completely
-				sky_bottom_override,   // double sky_temp_override,     // really cold average seed - ignore lowest_sky_row filter
-				sky_override_shrink,   // int    shrink_for_temp,       // shrink before finding hottest sky
-				sky_highest_min,       //  =   100; // lowest absolute value should not be higher (requires photometric)
-				clouds_en,             // enable clouds in the sky detection / processing
-				clouds_fom,            // maximal FOM for clouds (must be <=)         
-				clouds_spread,         // maximal spread for clouds (must be <=)         
-				clouds_disparity,      // maximal disparity for strong clouds         
-                clouds_weak,           // maximal strength for near definite clouds         		
-                clouds_strength,       // minimal strength for far strong clouds (definitely cloud)         		
-                clouds_not_strength,   // maximal strength for near maybe clouds (if it has strong cloud neighbor)         		
-				clouds_strong,         // true; // allow weak cloud if it has strong (1.5x) cloud neib
-				
-				wall_en,               // enable smooth walls detection/processing
-				wall_dflt,             // default (undetected) is wall (false - sky)
-				wall_str,              // minimal strength of the far object (small - just non-NaN disparity)         
-				wall_far,              // maximal disparity to consider cluster to be sky         
-				wall_near,             // minimal disparity to consider cluster to be wall         
-				
-			    treeline_en,           // enable treeline processing
-			    treeline_wide,         // look not only under, but diagonal too.
-			    treeline_height,       // maximal height of the treeline (tiles)
-			    treeline_width,        // minimal horizontal width of the treeline (tiles)
-			    treeline_lim_high,     // limit too high treeline (false - delete completely)
-			    treeline_str,          // treeline minimal strength 
-			    treeline_far,          // treeline min disparity (pix)
-			    treeline_near,         // treeline max disparity (pix)
-			    treeline_fg_str,       // pre-treeline FG objects (such as flat ground) minimal strength 
-			    treeline_fg_far,       // pre-treeline FG objects  min disparity (pix)
-			    treeline_fg_near,      // pre-treeline FG objects  max disparity (pix)
-				
-				indoors_en,            // true; // allow weak cloud if it has strong (1.5x) cloud neib
-				indoors_str,           //  0.5; // minimal strength of the far object         
-				indoors_disp,          //  0.8; // maximal minimal outdoor strong disparity         
-				indoors_min_out,       //  10;   // minimal strong far tiles to deny indoors
-				disp_boost_min,        // double disp_boost_min,    //  = 0.5;
-				disp_boost_diff,       //double disp_boost_diff,   //  = 0.35;
-				disp_boost_neibs,   //int    disp_boost_neibs,  //  = 2;
-				disp_boost_amount,  //double disp_boost_amount, //  = 2.0;
-				dsi[TwoQuadCLT.DSI_STRENGTH_AUX], // double [] strength,
-				dsi[TwoQuadCLT.DSI_SPREAD_AUX], // double [] spread,
-				dsi[TwoQuadCLT.DSI_DISPARITY_AUX_LMA], // double [] spread,
-				dsi[TwoQuadCLT.DSI_AVGVAL_AUX],//	double [] avg_val,
-				dbg_scene,          // QuadCLT   dbg_scene,    // use to save debug images if not null
-				batch_mode? -1: 1); /// debugLevel);        // int debugLevel)
+		
+		if (sky_extract) {
+			quadCLT_ref.setBlueSky  ( // initial BS from single scene
+					max_disparity,
+					max_disparity_strength,
+					sky_seed,           // double sky_seed, //  =       7.0;  // start with product of strength by diff_second below this
+					lma_seed,
+					seed_temp,     //double seed_temp, //         0.5;  // seed colder that this point between min and max temp						
+					sky_lim,            // double sky_lim, //   =      15.0; // then expand to product of strength by diff_second below this
+					lim_temp, //						double lim_temp, //          0.5;  // sky colder that this point between min and max temp						
+					sky_shrink,         // int    sky_shrink, //  =       4;
+					sky_expand_extra,   // int    sky_expand_extra, //  = 100; // 1?
+					sky_bottleneck,     //int    sky_bottleneck, // 
+					sky_reexpand_extra,   //int    sky_reexpand_extra,  // 9; re-expand after bottleneck in addition to how it was shrank
+					cold_scale, // =       0.2;  // <=1.0. 1.0 - disables temperature dependence
+					cold_frac, // =        0.005; // this and lower will scale fom by  cold_scale
+					hot_frac, // =         0.9;    // this and above will scale fom by 1.0
+					min_strength, // =     0.08;
+					seed_rows, // =        5; // sky should appear in this top rows 
+					lowest_sky_row,        //  =   50;// appears that low - invalid, remove completely
+					sky_bottom_override,   // double sky_temp_override,     // really cold average seed - ignore lowest_sky_row filter
+					sky_override_shrink,   // int    shrink_for_temp,       // shrink before finding hottest sky
+					sky_highest_min,       //  =   100; // lowest absolute value should not be higher (requires photometric)
+					clouds_en,             // enable clouds in the sky detection / processing
+					clouds_fom,            // maximal FOM for clouds (must be <=)         
+					clouds_spread,         // maximal spread for clouds (must be <=)         
+					clouds_disparity,      // maximal disparity for strong clouds         
+					clouds_weak,           // maximal strength for near definite clouds         		
+					clouds_strength,       // minimal strength for far strong clouds (definitely cloud)         		
+					clouds_not_strength,   // maximal strength for near maybe clouds (if it has strong cloud neighbor)         		
+					clouds_strong,         // true; // allow weak cloud if it has strong (1.5x) cloud neib
+
+					wall_en,               // enable smooth walls detection/processing
+					wall_dflt,             // default (undetected) is wall (false - sky)
+					wall_str,              // minimal strength of the far object (small - just non-NaN disparity)         
+					wall_far,              // maximal disparity to consider cluster to be sky         
+					wall_near,             // minimal disparity to consider cluster to be wall         
+
+					treeline_en,           // enable treeline processing
+					treeline_wide,         // look not only under, but diagonal too.
+					treeline_height,       // maximal height of the treeline (tiles)
+					treeline_width,        // minimal horizontal width of the treeline (tiles)
+					treeline_lim_high,     // limit too high treeline (false - delete completely)
+					treeline_str,          // treeline minimal strength 
+					treeline_far,          // treeline min disparity (pix)
+					treeline_near,         // treeline max disparity (pix)
+					treeline_fg_str,       // pre-treeline FG objects (such as flat ground) minimal strength 
+					treeline_fg_far,       // pre-treeline FG objects  min disparity (pix)
+					treeline_fg_near,      // pre-treeline FG objects  max disparity (pix)
+
+					indoors_en,            // true; // allow weak cloud if it has strong (1.5x) cloud neib
+					indoors_str,           //  0.5; // minimal strength of the far object         
+					indoors_disp,          //  0.8; // maximal minimal outdoor strong disparity         
+					indoors_min_out,       //  10;   // minimal strong far tiles to deny indoors
+					disp_boost_min,        // double disp_boost_min,    //  = 0.5;
+					disp_boost_diff,       //double disp_boost_diff,   //  = 0.35;
+					disp_boost_neibs,   //int    disp_boost_neibs,  //  = 2;
+					disp_boost_amount,  //double disp_boost_amount, //  = 2.0;
+					dsi[TwoQuadCLT.DSI_STRENGTH_AUX], // double [] strength,
+					dsi[TwoQuadCLT.DSI_SPREAD_AUX], // double [] spread,
+					dsi[TwoQuadCLT.DSI_DISPARITY_AUX_LMA], // double [] spread,
+					dsi[TwoQuadCLT.DSI_AVGVAL_AUX],//	double [] avg_val,
+					dbg_scene,          // QuadCLT   dbg_scene,    // use to save debug images if not null
+					batch_mode? -1: 1); /// debugLevel);        // int debugLevel)
+		} else {
+			quadCLT_ref.setBlueSky (new double [quadCLT_ref.tp.getTilesX()*quadCLT_ref.tp.getTilesY()]);
+		}
 		if (ran_photo_each) {
 //			quadCLT_ref.setBlueSky(null); // Reset blue sky - is it needed?
 			// see if blue sky was detected - rerun photoEach
@@ -4178,7 +4206,7 @@ public class OpticalFlow {
 						updateStatus,        // final boolean                                 updateStatus,
 						debugLevel);         // final int                                     debugLevel)			
 
-			} else {
+			} else if (sky_extract){
 				System.out.println("**** FAILED to detect non-empty Blue Sky after initial DSI in "+quadCLT_ref.getImageName()+
 						", if there is actual sky in the field of view, some \"sky\" parameters need to be adjusted.");
 				
@@ -4264,6 +4292,7 @@ public class OpticalFlow {
 	    double   treeline_fg_str =   clt_parameters.imp.treeline_fg_str;  //  0.8;  // pre-treeline FG objects (such as flat ground) minimal strength 
 	    double   treeline_fg_far =   clt_parameters.imp.treeline_fg_far;  //  0.2;  // pre-treeline FG objects  min disparity (pix)
 	    double   treeline_fg_near =  clt_parameters.imp.treeline_fg_near; //  0.5;  // pre-treeline FG objects  max disparity (pix)
+    	boolean sky_extract =        clt_parameters.imp.sky_extract;
 
 	    //	suspecting indoors (disabling sky)
 		boolean indoors_en =         clt_parameters.imp.indoors_en;       // true; // allow weak cloud if it has strong (1.5x) cloud neib
@@ -4297,7 +4326,6 @@ public class OpticalFlow {
 		if (dsi[TwoQuadCLT.DSI_SPREAD_AUX] == null) {
 			System.out.println("DSI_MAIN file has old format and does not have spread data, will recalculate.");
 		} else {
-//			boolean used_combo = false;
 			if ((dsi[TwoQuadCLT.DSI_BLUE_SKY_AUX] == null) || sky_recalc) { //
 				// Sets quadCLT_ref.dsi and blue sky (if exists)
 				combo_dsi = quadCLT_ref.restoreComboDSI(true); // result is full length, missing slices are null
@@ -4305,98 +4333,97 @@ public class OpticalFlow {
 					dsi[TwoQuadCLT.DSI_STRENGTH_AUX] =      combo_dsi[COMBO_DSN_INDX_STRENGTH];     // double [] strength,
 					dsi[TwoQuadCLT.DSI_DISPARITY_AUX_LMA]=  combo_dsi[COMBO_DSN_INDX_LMA]; //double [] disp_lma,
 					dsi[TwoQuadCLT.DSI_BLUE_SKY_AUX] =      combo_dsi[COMBO_DSN_INDX_BLUE_SKY]; // if exist, already set to quadCLT_ref.dsi 
-//					used_combo = true;
-				}
-				/*
-				if (quadCLT_ref.restoreInterDSI(true) > 0) {
-					double [][] combo_dsi = quadCLT_ref.dsi;
-					dsi[TwoQuadCLT.DSI_STRENGTH_AUX] =      combo_dsi[COMBO_DSN_INDX_STRENGTH];     // double [] strength,
-					dsi[TwoQuadCLT.DSI_DISPARITY_AUX_LMA]=  combo_dsi[COMBO_DSN_INDX_LMA]; //double [] disp_lma,
-					dsi[TwoQuadCLT.DSI_BLUE_SKY_AUX] =      combo_dsi[COMBO_DSN_INDX_BLUE_SKY];
-					used_combo = true;
 				}
-				*/
 			}
 			double [] ref_blue_sky = dsi[TwoQuadCLT.DSI_BLUE_SKY_AUX];
 			if ((ref_blue_sky == null) || sky_recalc) {
 				if (debugLevel > -3) {
 					System.out.println("Blue Sky does not exist or recalculation is forced. Calculating and updating photometrics");
 				}
-				double [] bs_strength = dsi[TwoQuadCLT.DSI_STRENGTH_AUX].clone();
-				// try getting combo dsi
-				if ((combo_dsi != null) && (scale_combo_strength != 1.0) && dsi[TwoQuadCLT.DSI_STRENGTH_AUX]!= null) {
-					for (int i = 0; i < dsi[TwoQuadCLT.DSI_STRENGTH_AUX].length; i++) {
-						bs_strength[i] *= scale_combo_strength;
-					}
-					//scale_combo_strength
-				}
-				quadCLT_ref.setBlueSky  (
-						max_disparity,
-						max_disparity_strength,
-						sky_seed,                              // double sky_seed, //  =       7.0;  // start with product of strength by diff_second below this
-						lma_seed,                              //          2.0;  // seed - disparity_lma limit
-						seed_temp,     //double seed_temp, //         0.5;  // seed colder that this point between min and max temp						
-						sky_lim,                               // double sky_lim, //   =      15.0; // then expand to product of strength by diff_second below this
-						lim_temp, //						double lim_temp, //          0.5;  // sky colder that this point between min and max temp						
-						sky_shrink,                            // int    sky_shrink, //  =       4;
-						sky_expand_extra,                      // int    sky_expand_extra, //  = 100; // 1?
-						sky_bottleneck,       //int    sky_bottleneck, // 
-						sky_reexpand_extra,   //int    sky_reexpand_extra,  // 9; re-expand after bottleneck in addition to how it was shrank
-						cold_scale, // =       0.2;  // <=1.0. 1.0 - disables temperature dependence
-						cold_frac, // =        0.005; // this and lower will scale fom by  cold_scale
-						hot_frac, // =         0.9;    // this and above will scale fom by 1.0
-						min_strength, // =     0.08;
-						seed_rows, // =        5; // sky should appear in this top rows
-						lowest_sky_row,        //  =     50;// appears that low - invalid, remove completely
-						sky_bottom_override,   // double sky_temp_override,     // really cold average seed - ignore lowest_sky_row filter
-						sky_override_shrink,   // int    shrink_for_temp,       // shrink before finding hottest sky
-						sky_highest_min,       //  =   100; // lowest absolute value should not be higher (requires photometric)
-
-						clouds_en,             // enable clouds in the sky detection / processing
-						clouds_fom,            // maximal FOM for clouds (must be <=)         
-						clouds_spread,         // maximal spread for clouds (must be <=)         
-						clouds_disparity,      // maximal disparity for strong clouds         
-		                clouds_weak,           // maximal strength for near definite clouds         		
-		                clouds_strength,       // minimal strength for far strong clouds (definitely cloud)         		
-		                clouds_not_strength,   // maximal strength for near maybe clouds (if it has strong cloud neighbor)         		
-						clouds_strong,         // true; // allow weak cloud if it has strong (1.5x) cloud neib
-
-						wall_en,               // enable smooth walls detection/processing
-						wall_dflt,             // default (undetected) is wall (false - sky)
-						wall_str,              // minimal strength of the far object (small - just non-NaN disparity)         
-						wall_far,              // maximal disparity to consider cluster to be sky         
-						wall_near,             // minimal disparity to consider cluster to be wall         
-
-					    treeline_en,           // enable treeline processing
-					    treeline_wide,         // look not only under, but diagonal too.
-					    treeline_height,       // maximal height of the treeline (tiles)
-					    treeline_width,        // minimal horizontal width of the treeline (tiles)
-					    treeline_lim_high,     // limit too high treeline (false - delete completely)
-					    treeline_str,          // treeline minimal strength 
-					    treeline_far,          // treeline min disparity (pix)
-					    treeline_near,         // treeline max disparity (pix)
-					    treeline_fg_str,       // pre-treeline FG objects (such as flat ground) minimal strength 
-					    treeline_fg_far,       // pre-treeline FG objects  min disparity (pix)
-					    treeline_fg_near,      // pre-treeline FG objects  max disparity (pix)
-
-					    indoors_en,            // true; // allow weak cloud if it has strong (1.5x) cloud neib
-						indoors_str,           //  0.5; // minimal strength of the far object         
-						indoors_disp,          //  0.8; // maximal minimal outdoor strong disparity         
-						indoors_min_out,       //  10;   // minimal strong far tiles to deny indoors
-						
-						disp_boost_min,        // double disp_boost_min,    //  = 0.5;
-						disp_boost_diff,       //double disp_boost_diff,   //  = 0.35;
-						disp_boost_neibs,      //int    disp_boost_neibs,  //  = 2;
-						disp_boost_amount,     //double disp_boost_amount, //  = 2.0;
-						bs_strength, // dsi[TwoQuadCLT.DSI_STRENGTH_AUX],      // double [] strength,
-						dsi[TwoQuadCLT.DSI_SPREAD_AUX],        // double [] spread,
-						dsi[TwoQuadCLT.DSI_DISPARITY_AUX_LMA], //double [] disp_lma,
-						dsi[TwoQuadCLT.DSI_AVGVAL_AUX],//	double [] avg_val,
-						dbg_scene,          // QuadCLT   dbg_scene,    // use to save debug images if not null
-						debugLevel);        // int debugLevel)
-				if (debugLevel > -3) {
-					System.out.println("Calculated missing Blue Sky in "+quadCLT_ref.getImageName()+
-							", re-running photoEach()");
+				if (sky_extract) { 
+					double [] bs_strength = dsi[TwoQuadCLT.DSI_STRENGTH_AUX].clone();
+					// try getting combo dsi
+					if ((combo_dsi != null) && (scale_combo_strength != 1.0) && dsi[TwoQuadCLT.DSI_STRENGTH_AUX]!= null) {
+						for (int i = 0; i < dsi[TwoQuadCLT.DSI_STRENGTH_AUX].length; i++) {
+							bs_strength[i] *= scale_combo_strength;
+						}
+						//scale_combo_strength
+					}
+					quadCLT_ref.setBlueSky  (
+							max_disparity,
+							max_disparity_strength,
+							sky_seed,                              // double sky_seed, //  =       7.0;  // start with product of strength by diff_second below this
+							lma_seed,                              //          2.0;  // seed - disparity_lma limit
+							seed_temp,     //double seed_temp, //         0.5;  // seed colder that this point between min and max temp						
+							sky_lim,                               // double sky_lim, //   =      15.0; // then expand to product of strength by diff_second below this
+							lim_temp, //						double lim_temp, //          0.5;  // sky colder that this point between min and max temp						
+							sky_shrink,                            // int    sky_shrink, //  =       4;
+							sky_expand_extra,                      // int    sky_expand_extra, //  = 100; // 1?
+							sky_bottleneck,       //int    sky_bottleneck, // 
+							sky_reexpand_extra,   //int    sky_reexpand_extra,  // 9; re-expand after bottleneck in addition to how it was shrank
+							cold_scale, // =       0.2;  // <=1.0. 1.0 - disables temperature dependence
+							cold_frac, // =        0.005; // this and lower will scale fom by  cold_scale
+							hot_frac, // =         0.9;    // this and above will scale fom by 1.0
+							min_strength, // =     0.08;
+							seed_rows, // =        5; // sky should appear in this top rows
+							lowest_sky_row,        //  =     50;// appears that low - invalid, remove completely
+							sky_bottom_override,   // double sky_temp_override,     // really cold average seed - ignore lowest_sky_row filter
+							sky_override_shrink,   // int    shrink_for_temp,       // shrink before finding hottest sky
+							sky_highest_min,       //  =   100; // lowest absolute value should not be higher (requires photometric)
+
+							clouds_en,             // enable clouds in the sky detection / processing
+							clouds_fom,            // maximal FOM for clouds (must be <=)         
+							clouds_spread,         // maximal spread for clouds (must be <=)         
+							clouds_disparity,      // maximal disparity for strong clouds         
+							clouds_weak,           // maximal strength for near definite clouds         		
+							clouds_strength,       // minimal strength for far strong clouds (definitely cloud)         		
+							clouds_not_strength,   // maximal strength for near maybe clouds (if it has strong cloud neighbor)         		
+							clouds_strong,         // true; // allow weak cloud if it has strong (1.5x) cloud neib
+
+							wall_en,               // enable smooth walls detection/processing
+							wall_dflt,             // default (undetected) is wall (false - sky)
+							wall_str,              // minimal strength of the far object (small - just non-NaN disparity)         
+							wall_far,              // maximal disparity to consider cluster to be sky         
+							wall_near,             // minimal disparity to consider cluster to be wall         
+
+							treeline_en,           // enable treeline processing
+							treeline_wide,         // look not only under, but diagonal too.
+							treeline_height,       // maximal height of the treeline (tiles)
+							treeline_width,        // minimal horizontal width of the treeline (tiles)
+							treeline_lim_high,     // limit too high treeline (false - delete completely)
+							treeline_str,          // treeline minimal strength 
+							treeline_far,          // treeline min disparity (pix)
+							treeline_near,         // treeline max disparity (pix)
+							treeline_fg_str,       // pre-treeline FG objects (such as flat ground) minimal strength 
+							treeline_fg_far,       // pre-treeline FG objects  min disparity (pix)
+							treeline_fg_near,      // pre-treeline FG objects  max disparity (pix)
+
+							indoors_en,            // true; // allow weak cloud if it has strong (1.5x) cloud neib
+							indoors_str,           //  0.5; // minimal strength of the far object         
+							indoors_disp,          //  0.8; // maximal minimal outdoor strong disparity         
+							indoors_min_out,       //  10;   // minimal strong far tiles to deny indoors
+
+							disp_boost_min,        // double disp_boost_min,    //  = 0.5;
+							disp_boost_diff,       //double disp_boost_diff,   //  = 0.35;
+							disp_boost_neibs,      //int    disp_boost_neibs,  //  = 2;
+							disp_boost_amount,     //double disp_boost_amount, //  = 2.0;
+							bs_strength, // dsi[TwoQuadCLT.DSI_STRENGTH_AUX],      // double [] strength,
+							dsi[TwoQuadCLT.DSI_SPREAD_AUX],        // double [] spread,
+							dsi[TwoQuadCLT.DSI_DISPARITY_AUX_LMA], //double [] disp_lma,
+							dsi[TwoQuadCLT.DSI_AVGVAL_AUX],//	double [] avg_val,
+							dbg_scene,          // QuadCLT   dbg_scene,    // use to save debug images if not null
+							debugLevel);        // int debugLevel)
+					if (debugLevel > -3) {
+						System.out.println("Calculated missing Blue Sky in "+quadCLT_ref.getImageName()+
+								", re-running photoEach()");
+					}
+				} else {
+					quadCLT_ref.setBlueSky (new double [quadCLT_ref.tp.getTilesX()*quadCLT_ref.tp.getTilesY()]);
+					if (debugLevel > -3) {
+						System.out.println("Used zeros for the missing Blue Sky in "+quadCLT_ref.getImageName()+
+								", not re-running photoEach()");
+					}
+
 				}
 				// save inter-lma if available
 				if (combo_dsi != null) {
@@ -4414,19 +4441,20 @@ public class OpticalFlow {
 							quadCLT_ref.tp.getTilesX(),   // tilesX,                // int         width,
 							quadCLT_ref.tp.getTilesY());  // int         height)
 				}
-
-				photoEach(
-						clt_parameters,      // CLTParameters                                 clt_parameters,
-						colorProcParameters, // ColorProcParameters                           colorProcParameters,
-						quadCLT_main,        // QuadCLT                                       quadCLT_main, // tiles should be set
-						quadCLT_ref,         // QuadCLT                                       quadCLT_ref, // tiles should be set
-						dsi,                 // final double [][]                             dsi,
-						// just once?
-						clt_parameters.photo_num_full, // int      photo_num_full =              
-						batch_mode,          // final boolean                                 batch_mode,
-						threadsMax,          // final int                                     threadsMax,  // int               threadsMax,
-						updateStatus,        // final boolean                                 updateStatus,
-						debugLevel);         // final int                                     debugLevel)			
+				if (sky_extract) { 
+					photoEach(
+							clt_parameters,      // CLTParameters                                 clt_parameters,
+							colorProcParameters, // ColorProcParameters                           colorProcParameters,
+							quadCLT_main,        // QuadCLT                                       quadCLT_main, // tiles should be set
+							quadCLT_ref,         // QuadCLT                                       quadCLT_ref, // tiles should be set
+							dsi,                 // final double [][]                             dsi,
+							// just once?
+							clt_parameters.photo_num_full, // int      photo_num_full =              
+							batch_mode,          // final boolean                                 batch_mode,
+							threadsMax,          // final int                                     threadsMax,  // int               threadsMax,
+							updateStatus,        // final boolean                                 updateStatus,
+							debugLevel);         // final int                                     debugLevel)	
+				}
 			} else {
 				if (debugLevel > -3) {
 					System.out.println("Blue Sky is available, reusing it, no need for updating photometrics");
@@ -4442,6 +4470,7 @@ public class OpticalFlow {
     		String              set_name,
 			QuadCLT[]           quadCLTs,
 			int                 earliest_scene,
+			int                 last_scene,
 			int                 ref_index,
 			double [][]         combo_dsn_final_filtered,
 			int                 threadsMax,  // maximal number of threads to launch
@@ -4494,7 +4523,7 @@ public class OpticalFlow {
 ///			ErsCorrection ers_reference = quadCLTs[ref_index].getErsCorrection();
 		}
 		// propagate to other scenes in this sequence
-		for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
+		for (int scene_index =  last_scene; scene_index >= earliest_scene ; scene_index--) if (scene_index != ref_index){
 			//						quadCLTs[scene_index] = (QuadCLT) quadCLT_main.spawnNoModelQuadCLT( // restores image data
 			// to include ref scene photometric calibration
 			/*
@@ -4519,9 +4548,10 @@ public class OpticalFlow {
 	
 	
 	/**
-	 * Build series of poses from just a single (reference) scene.
+	 * 
+	 * @param batch_mode
 	 * @param quadCLT_main
-	 * @param ref_index
+	 * @param last_index
 	 * @param clt_parameters
 	 * @param debayerParameters
 	 * @param colorProcParameters
@@ -4530,16 +4560,19 @@ public class OpticalFlow {
 	 * @param equirectangularParameters
 	 * @param properties
 	 * @param reset_from_extrinsics
+	 * @param videos
+	 * @param stereo_widths
+	 * @param start_ref_pointers sho
 	 * @param threadsMax
 	 * @param updateStatus
 	 * @param debugLevel
-	 * @return null on failure or after LY adjustments, or String path to the reference model directory 
+	 * @return
 	 * @throws Exception
 	 */
     public String buildSeries(
     		boolean                                              batch_mode,
 			QuadCLT                                              quadCLT_main, // tiles should be set
-			int                                                  ref_index, // -1 - last
+			int                                                  last_index, // -1 - last
 			CLTParameters                                        clt_parameters,
 			EyesisCorrectionParameters.DebayerParameters         debayerParameters,
 			ColorProcParameters                                  colorProcParameters,
@@ -4557,7 +4590,8 @@ public class OpticalFlow {
 			final int        debugLevel)  throws Exception
     {
     	int earliest_scene = 0;      // increase on failure
-    	
+    	boolean center_reference =           clt_parameters.imp.center_reference;
+    	int min_num_scenes =                 clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
     	boolean build_ref_dsi =              clt_parameters.imp.force_ref_dsi;
     	boolean force_initial_orientations = clt_parameters.imp.force_orientations ;
     	boolean run_ly =                     clt_parameters.imp.run_ly; // run LY adjust and exit
@@ -4574,7 +4608,7 @@ public class OpticalFlow {
     	
     	boolean show_dsi_image =             clt_parameters.imp.show_ranges && !batch_mode;
     	boolean show_images =                clt_parameters.imp.show_images && !batch_mode;
-    	int min_num_scenes =                 clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
+//    	int min_num_scenes =                 clt_parameters.imp.min_num_scenes; // abandon series if there are less than this number of scenes in it 
     	int max_num_scenes =                 clt_parameters.imp.max_num_scenes; // cut longer series 
     	
     	boolean show_images_bgfg =           clt_parameters.imp.show_images_bgfg && !batch_mode;
@@ -4681,74 +4715,91 @@ public class OpticalFlow {
 		}
 		// set_channels will include all 99 scenes even as quadCLTs.length matches ref_index
 		QuadCLT.SetChannels [] set_channels=quadCLT_main.setChannels(debugLevel); 
-		if (ref_index < 0) {
-			ref_index += set_channels.length; 
+		if (last_index < 0) {
+			last_index += set_channels.length; 
 		}
 		if (start_ref_pointers != null) {
 			start_ref_pointers[0] = 0;
-			start_ref_pointers[1] = ref_index;
+			start_ref_pointers[1] = last_index;
 		}
 		
-		QuadCLT [] quadCLTs = new QuadCLT [ref_index+1]; //  [set_channels.length];
+		QuadCLT [] quadCLTs = new QuadCLT [last_index+1]; //  [set_channels.length];
 		
 		//start_index
 		double [][][] scenes_xyzatr =      new double [quadCLTs.length][][]; // previous scene relative to the next one
-		scenes_xyzatr[ref_index] = new double[2][3]; // all zeros
+		scenes_xyzatr[last_index] = new double[2][3]; // all zeros
 		// See if build_ref_dsi is needed
 		if (!build_ref_dsi) {
-			quadCLTs[ref_index] = (QuadCLT) quadCLT_main.spawnNoModelQuadCLT( // will conditionImageSet
-					set_channels[ref_index].set_name,
+			quadCLTs[last_index] = (QuadCLT) quadCLT_main.spawnNoModelQuadCLT( // will conditionImageSet
+					set_channels[last_index].set_name,
 					clt_parameters,
 					colorProcParameters, //
 					threadsMax,
 					debugLevel-2);
-			if ((quadCLTs[ref_index] == null) || !quadCLTs[ref_index].dsiExists()) {
-				if (debugLevel >-2) {
-					System.out.println("DSI data for scene "+set_channels[ref_index].set_name+" does not exist, forcing its calculation.");
+			if (center_reference && (quadCLTs[last_index] != null)) {
+				quadCLTs[last_index].restoreInterProperties( // restore properties for interscene processing (extrinsics, ers, ...) // get relative poses (98)
+						null, // String path,             // full name with extension or null to use x3d directory
+						false, // boolean all_properties,//				null, // Properties properties,   // if null - will only save extrinsics)
+						debugLevel);
+				int center_index =quadCLTs[last_index].getReferenceIndex(null);
+				if (center_index == -1) {
+					build_ref_dsi = true;
+				}
+			} else {
+				if ((quadCLTs[last_index] == null) || !quadCLTs[last_index].dsiExists()) {
+					if (debugLevel >-2) {
+						System.out.println("DSI data for scene "+set_channels[last_index].set_name+" does not exist, forcing its calculation.");
+					}
+					build_ref_dsi = true;
 				}
-				build_ref_dsi = true;
 			}
 		}
-		if (!build_ref_dsi && (quadCLTs[ref_index] != null)) {
-			quadCLTs[ref_index].restoreInterProperties(null, false, debugLevel); //null
+		if (!build_ref_dsi && (quadCLTs[last_index] != null)) {
+			quadCLTs[last_index].restoreInterProperties(null, false, debugLevel); //null
 		}
 		// 1. Reference scene DSI
 //		while ((quadCLTs[ref_index] == null) || (quadCLTs[ref_index].getBlueSky() == null)) { // null
-		while ((quadCLTs[ref_index] == null) || !quadCLTs[ref_index].hasBlueSky()) { // null
+//		String center_ts = null; // name of the center scene - not yet known
+		while ((quadCLTs[last_index] == null) || !quadCLTs[last_index].hasBlueSky()) { // null
 			//hasBlueSky()
 			if (build_ref_dsi) {
 				TwoQuadCLT.copyJP4src( // actually there is no sense to process multiple image sets. Combine with other
 						// processing?
-						set_channels[ref_index].set_name, // String set_name
+						set_channels[last_index].set_name, // String set_name
 						quadCLT_main, // QuadCLT quadCLT_main,
 						null, // QuadCLT quadCLT_aux,
 						clt_parameters, // EyesisCorrectionParameters.DCTParameters dct_parameters,
 						true, // boolean                                  skip_existing,
 						true, // false, // boolean                                  search_KML,
 						debugLevel);
-				quadCLTs[ref_index] = (QuadCLT) quadCLT_main.spawnNoModelQuadCLT( // will conditionImageSet
-						set_channels[ref_index].set_name,
+				quadCLTs[last_index] = (QuadCLT) quadCLT_main.spawnNoModelQuadCLT( // will conditionImageSet
+						set_channels[last_index].set_name,
 						clt_parameters,
 						colorProcParameters, //
 						threadsMax,
 						debugLevel-2);
-				quadCLTs[ref_index].saveQuadClt(); // to re-load new set of Bayer images to the GPU (do nothing for CPU) and Geometry
+				quadCLTs[last_index].saveQuadClt(); // to re-load new set of Bayer images to the GPU (do nothing for CPU) and Geometry
+				boolean fast = center_reference; // false;
+				if (center_reference) { // invalidate pointer to recalculate
+					quadCLTs[last_index].setRefPointer((String) null); 
+				}
 //				quadCLTs[ref_index] = 
 				buildRefDSI( // returned is a different instance than input -FIXED
 						clt_parameters,                   // CLTParameters                                 clt_parameters,
+						fast,                             // boolean  fast,
 						debayerParameters,                // EyesisCorrectionParameters.DebayerParameters  debayerParameters,
 						colorProcParameters,              // ColorProcParameters                           colorProcParameters,
 					    rgbParameters, // EyesisCorrectionParameters.RGBParameters      rgbParameters,
 						batch_mode,                       // boolean                             batch_mode,
-						set_channels[ref_index].set_name, // String                              set_name,
+						set_channels[last_index].set_name, // String                              set_name,
 						quadCLT_main,                     // QuadCLT                             quadCLT_main, // tiles should be set
-						quadCLTs[ref_index],              // QuadCLT                             quadCLT_ref, // tiles should be set
+						quadCLTs[last_index],              // QuadCLT                             quadCLT_ref, // tiles should be set
 					    threadsMax,                       // final int                                     threadsMax,  // maximal number of threads to launch
 					    updateStatus,                     // final boolean    updateStatus,
 					    debugLevel);                      // final int        debugLevel);
 				earliest_scene = 0;  // reset failures, try to use again all scenes
 				// Copy source files for other scenes
-				for (int scene_index =  ref_index - 1; scene_index >= 0 ; scene_index--) {
+				for (int scene_index =  last_index - 1; scene_index >= 0 ; scene_index--) {
 					TwoQuadCLT.copyJP4src( // actually there is no sense to process multiple image sets. Combine with other
 							// processing?
 							set_channels[scene_index].set_name, // String set_name
@@ -4759,119 +4810,199 @@ public class OpticalFlow {
 							true, // false, // boolean                                  search_KML,
 							debugLevel);
 				} // split cycles to remove output clutter
-				
-				
-				
 			} else {// if (build_ref_dsi) {
 				// need to read photometric from reference scene -INTERFRAME.corr-xml, and if it exists - set and propagate to main?
+				// try to run reuseRefDSI() even for center? 
 				reuseRefDSI(
 						clt_parameters,      // CLTParameters  clt_parameters,
 						colorProcParameters, // ColorProcParameters colorProcParameters,
 						quadCLT_main,        // QuadCLT        quadCLT_main, // tiles should be set
-						quadCLTs[ref_index], // QuadCLT        quadCLT_ref, // tiles should be set
+						quadCLTs[last_index], // QuadCLT        quadCLT_ref, // tiles should be set
 						batch_mode,			 // final boolean       batch_mode,
 						threadsMax, 		 // final int           threadsMax,
 						updateStatus,		 // final boolean       updateStatus,
 						debugLevel);         // int            debugLevel)
+				if (center_reference) {
+//					center_ts = getReferenceTimestamp();
+				}
 			}
-			
 		} // while (blue_sky == null) 
 		
-//		ref_blue_sky = quadCLTs[ref_index].getBlueSky();
-		ref_blue_sky = quadCLTs[ref_index].getDoubleBlueSky();
+		// by now if it was center_reference the last ref_dsi is restored with optional blue sky.
 		
-		quadCLTs[ref_index] = (QuadCLT) quadCLT_main.spawnQuadCLT( // restores dsi from "DSI-MAIN"
-				set_channels[ref_index].set_name,
+		
+		ref_blue_sky = quadCLTs[last_index].getDoubleBlueSky();
+		quadCLTs[last_index] = (QuadCLT) quadCLT_main.spawnQuadCLT( // restores dsi from "DSI-MAIN"
+				set_channels[last_index].set_name,
 				clt_parameters,
 				colorProcParameters, //
 				threadsMax,
 				debugLevel);
-		quadCLTs[ref_index].setQuadClt(); // just in case ?
-		quadCLTs[ref_index].setBlueSky(ref_blue_sky); //quadCLTs[ref_index].dsi has it 
-		
-		quadCLTs[ref_index].setDSRBG(
+		quadCLTs[last_index].setQuadClt(); // just in case ?
+		quadCLTs[last_index].setBlueSky(ref_blue_sky); //quadCLTs[ref_index].dsi has it 
+		quadCLTs[last_index].setDSRBG(
 				clt_parameters, // CLTParameters  clt_parameters,
 				threadsMax,     // int            threadsMax,  // maximal number of threads to launch
 				updateStatus,   // boolean        updateStatus,
 				debugLevel);    // int            debugLevel)
-		ErsCorrection ers_reference = quadCLTs[ref_index].getErsCorrection();
 		
 		if (!force_initial_orientations) {
-			if ((quadCLTs[ref_index].getNumOrient() == 0) ||(!quadCLTs[ref_index].propertiesContainString (".scenes_"))) {
-				if (debugLevel >-3) {
-					System.out.println("Egomotion data for scene "+set_channels[ref_index].set_name+" does not exist, forcing its calculation.");
+			if (center_reference) {
+				int center_index =quadCLTs[last_index].getReferenceIndex(null);
+				if (center_index == -1) {
+					force_initial_orientations = true;
+				}
+
+			} else {
+				if ((quadCLTs[last_index].getNumOrient() == 0) ||(!quadCLTs[last_index].propertiesContainString (".scenes_"))) {
+					if (debugLevel >-3) {
+						System.out.println("Egomotion data for scene "+set_channels[last_index].set_name+" does not exist, forcing its calculation.");
+					}
+					force_initial_orientations = true;
 				}
-				force_initial_orientations = true;
 			}
 		}
 
 		// Build initial orientations
-		boolean ims_use =         clt_parameters.imp.ims_use;
+		boolean ims_use =      center_reference || clt_parameters.imp.ims_use; // center works with IMS only
+///		int center_index=ref_index;
+		int ref_index = last_index; // old versions
 		if (force_initial_orientations && !reuse_video) {
-			boolean OK =ims_use ?	
-					Interscene.setInitialOrientationsIms(
-							clt_parameters,      // final CLTParameters          clt_parameters,
-							colorProcParameters, // final ColorProcParameters    colorProcParameters,
-							quadCLTs,            //  final QuadCLT[]              quadCLTs, //
-							ref_index,           // final int                    ref_index,
-							set_channels,        // final QuadCLT.SetChannels [] set_channels,
-							batch_mode,          // final boolean                batch_mode,
-							earliest_scene,      // int                          earliest_scene,
-							start_ref_pointers,  // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
-							threadsMax,          // final int                    threadsMax,
-							updateStatus,        // final boolean                updateStatus,
-							debugLevel):         // final int                    debugLevel)
-								Interscene.setInitialOrientations(
-										clt_parameters,      // final CLTParameters          clt_parameters,
-										colorProcParameters, // final ColorProcParameters    colorProcParameters,
-										quadCLTs,            //  final QuadCLT[]              quadCLTs, //
-										ref_index,           // final int                    ref_index,
-										set_channels,        // final QuadCLT.SetChannels [] set_channels,
-										batch_mode,          // final boolean                batch_mode,
-										earliest_scene,      // int                          earliest_scene,
-										start_ref_pointers,  // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
-										threadsMax,          // final int                    threadsMax,
-										updateStatus,        // final boolean                updateStatus,
-										debugLevel);         // final int                    debugLevel)
+			boolean OK = false;
+			int es1 = -1;
+			if (center_reference) {
+				ref_index = Interscene.setInitialOrientationsCenterIms(
+						clt_parameters,      // final CLTParameters          clt_parameters,
+						min_num_scenes,      // int                          min_num_scenes,
+						colorProcParameters, // final ColorProcParameters    colorProcParameters,
+						debayerParameters,	 // EyesisCorrectionParameters.DebayerParameters debayerParameters,			
+						rgbParameters,       // EyesisCorrectionParameters.RGBParameters     rgbParameters,
+						quadCLT_main,        // QuadCLT                                      quadCLT_main, // tiles should be set			
+						quadCLTs,            // final QuadCLT[]              quadCLTs, //
+						last_index,           // final int                    ref_index,
+						set_channels,        // final QuadCLT.SetChannels [] set_channels,
+						batch_mode,          // final boolean                batch_mode,
+						earliest_scene,      // int                          earliest_scene,
+						start_ref_pointers,  // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
+						threadsMax,          // final int                    threadsMax,
+						updateStatus,        // final boolean                updateStatus,
+						debugLevel);
+				OK = ref_index >= 0;
+			} else if (ims_use) {	
+				es1 =	Interscene.setInitialOrientationsIms(
+						clt_parameters,      // final CLTParameters          clt_parameters,
+						min_num_scenes,      // int                          min_num_scenes,
+						colorProcParameters, // final ColorProcParameters    colorProcParameters,
+						quadCLTs,            //  final QuadCLT[]              quadCLTs, //
+						last_index,           // final int                    ref_index,
+						set_channels,        // final QuadCLT.SetChannels [] set_channels,
+						batch_mode,          // final boolean                batch_mode,
+						earliest_scene,      // int                          earliest_scene,
+						start_ref_pointers,  // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
+						threadsMax,          // final int                    threadsMax,
+						updateStatus,        // final boolean                updateStatus,
+						debugLevel);         // final int                    debugLevel)
+				OK = es1 >=0;
+			} else {
+				es1=	Interscene.setInitialOrientations(
+						clt_parameters,      // final CLTParameters          clt_parameters,
+						min_num_scenes,      // int                          min_num_scenes,
+						colorProcParameters, // final ColorProcParameters    colorProcParameters,
+						quadCLTs,            //  final QuadCLT[]              quadCLTs, //
+						last_index,           // final int                    ref_index,
+						set_channels,        // final QuadCLT.SetChannels [] set_channels,
+						batch_mode,          // final boolean                batch_mode,
+						earliest_scene,      // int                          earliest_scene,
+						start_ref_pointers,  // int []                       start_ref_pointers, // [0] - earliest valid scene, [1] ref_index
+						threadsMax,          // final int                    threadsMax,
+						updateStatus,        // final boolean                updateStatus,
+						debugLevel);         // final int                    debugLevel)
+				OK = es1 >=0;
+			}
 			if (!OK) {
 				return null;
 			}
 		} else {// if (build_orientations) {
 			if (!reuse_video) { // reuse_video only uses reference scene
-				for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
-					// to include ref scene photometric calibration
-					quadCLTs[scene_index] = quadCLTs[ref_index].spawnNoModelQuadCLT( // restores image data
-							set_channels[scene_index].set_name,
-							clt_parameters,
-							colorProcParameters, //
-							threadsMax,
-							debugLevel-2);
+//				if (ref_index)
+				for (int scene_index =  last_index; scene_index >= earliest_scene ; scene_index--) {
+					if (scene_index != ref_index) {
+						// to include ref scene photometric calibration
+						quadCLTs[scene_index] = quadCLTs[ref_index].spawnNoModelQuadCLT( // restores image data
+								set_channels[scene_index].set_name,
+								clt_parameters,
+								colorProcParameters, //
+								threadsMax,
+								debugLevel-2);
+					}
 				}
 			}
-			for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
-				if ((!quadCLTs[ref_index].tsExists(set_channels[scene_index].set_name)) ||
-						((ref_index - scene_index) >= max_num_scenes)) {
-					earliest_scene = scene_index + 1;
-				}
+			if (center_reference) {
+				ref_index =quadCLTs[last_index].getReferenceIndex(quadCLTs);
+				quadCLTs[ref_index].restoreInterProperties(null, false, debugLevel); //null
+				int [] first_last_index = quadCLTs[ref_index].getFirstLastIndex(quadCLTs);
+				earliest_scene = first_last_index[0];
+				last_index =  first_last_index[1];
 			}
-			if ((ref_index - earliest_scene + 1) < min_num_scenes) {
-				System.out.println("2.Total number of useful scenes = "+(ref_index - earliest_scene + 1)+
-						" < "+min_num_scenes+". Scrapping this series.");
-				if (start_ref_pointers != null) {
-					start_ref_pointers[0] = earliest_scene;
+			// TODO: 10.17.2023 - verify OK to remove this - next should be already tested during initial orientation
+			/*
+			if (!center_reference) { // for center reference it should be set during initial orientation. Or for all too? 
+				for (int scene_index =  ref_index - 1; scene_index >= earliest_scene ; scene_index--) {
+					if ((!quadCLTs[ref_index].tsExists(set_channels[scene_index].set_name)) ||
+							((ref_index - scene_index) >= max_num_scenes)) {
+						earliest_scene = scene_index + 1;
+					}
+				}
+				if ((ref_index - earliest_scene + 1) < min_num_scenes) {
+					System.out.println("2.Total number of useful scenes = "+(ref_index - earliest_scene + 1)+
+							" < "+min_num_scenes+". Scrapping this series.");
+					if (start_ref_pointers != null) {
+						start_ref_pointers[0] = earliest_scene;
+					}
+					return null;
 				}
-				return null;
+			}
+			*/
+		}
+		// both after initial orientations and if loaded
+		if (ref_index != last_index) { // do that for reference (center) scene again (was for the last scene)
+			// FIXME: Will deal with blue sky later - not needed for drones.
+///				ref_blue_sky = quadCLTs[ref_index].getDoubleBlueSky();
+			quadCLTs[ref_index] = (QuadCLT) quadCLT_main.spawnQuadCLT( // restores dsi from "DSI-MAIN"
+					set_channels[ref_index].set_name,
+					clt_parameters,
+					colorProcParameters, //
+					threadsMax,
+					debugLevel);
+			quadCLTs[ref_index].setQuadClt(); // just in case ?
+///				quadCLTs[ref_index].setBlueSky(ref_blue_sky); //quadCLTs[ref_index].dsi has it 
+			quadCLTs[ref_index].setDSRBG(
+					clt_parameters, // CLTParameters  clt_parameters,
+					threadsMax,     // int            threadsMax,  // maximal number of threads to launch
+					updateStatus,   // boolean        updateStatus,
+					debugLevel);    // int            debugLevel)
+			
+			int [] first_last_index = quadCLTs[ref_index].getFirstLastIndex(quadCLTs);
+			if (first_last_index == null) {
+				System.out.println("BUG: saved reference scene "+quadCLTs[ref_index].getImageName()+
+						" does not have first/last scenes defined");
+				System.out.println();
+			} else {
+				earliest_scene = first_last_index[0];
+				last_index =  first_last_index[1];
 			}
 		}
-		// just in case that orientations were calculated before:
-//		earliest_scene = getEarliestScene(quadCLTs);
 		
+		// just in case that orientations were calculated before:
+		//		earliest_scene = getEarliestScene(quadCLTs);
+
 		double [][] combo_dsn_final = null;
+		// below ref_index is not necessary the last (fix all where it is supposed to be the last
+		ErsCorrection ers_reference = quadCLTs[ref_index].getErsCorrection();
 		while (!reuse_video && ((quadCLTs[ref_index].getNumOrient() < min_num_orient) || (quadCLTs[ref_index].getNumAccum() < min_num_interscene))) {
 			if ((quadCLTs[ref_index].getNumAccum() < min_num_interscene) &&
 					((quadCLTs[ref_index].getNumAccum() <  quadCLTs[ref_index].getNumOrient())||
 							(quadCLTs[ref_index].getNumOrient() >= min_num_orient))) {
-				
 				boolean done_sfm = false;
 				if (quadCLTs[ref_index].getNumAccum() > 0) {
 					double mb_max_gain = clt_parameters.imp.mb_max_gain;
@@ -4889,34 +5020,43 @@ public class OpticalFlow {
 					*/
 //					int      num_avg_pairs = 16; // number of scene pairs to average
 					int      num_avg_pairs = 1;  // number of scene pairs to average
-					QuadCLT[][] scenes_pairs = new QuadCLT[num_avg_pairs][2];
+					QuadCLT[][][] scenes_seq_pairs = new QuadCLT[3][num_avg_pairs][2];
 					for (int i = 0; i < num_avg_pairs; i++) {
-//						scenes_pairs[i][0] = quadCLTs[ref_index - 1 - i];
-						scenes_pairs[i][0] = quadCLTs[ref_index - i];
-						scenes_pairs[i][1] = quadCLTs[earliest_scene + num_avg_pairs - 1 - i];
+						scenes_seq_pairs[0][i][0] = quadCLTs[ref_index + num_avg_pairs - 1 - i];
+						scenes_seq_pairs[0][i][1] = quadCLTs[earliest_scene + num_avg_pairs - 1 - i];
+						
+						scenes_seq_pairs[1][i][0] = quadCLTs[last_index - i];
+						scenes_seq_pairs[1][i][1] = quadCLTs[ref_index - i];
+						
+						scenes_seq_pairs[2][i][0] = quadCLTs[last_index - i];
+						scenes_seq_pairs[2][i][1] = quadCLTs[earliest_scene + num_avg_pairs - 1 - i];
 					}
-					double [][] sfm_dsn = StructureFromMotion.sfmPair(
-							clt_parameters,           // final CLTParameters          clt_parameters,
-							quadCLTs[ref_index],      // final QuadCLT                ref_scene,
-							scenes_pairs,             // final QuadCLT [][]           scenes_pairs,
-							mb_max_gain,              // double                       mb_max_gain,
-							batch_mode,               // final boolean                batch_mode,
-							debugLevel);              // final int                    debugLevel)
-					if (sfm_dsn != null) {
-						combo_dsn_final = sfm_dsn;
-						done_sfm = true;
+					
+					for (int num_seq = 0; num_seq < scenes_seq_pairs.length; num_seq++) {
+						double [][] sfm_dsn = StructureFromMotion.sfmPair(
+								clt_parameters,           // final CLTParameters          clt_parameters,
+								quadCLTs[ref_index],      // final QuadCLT                ref_scene,
+								scenes_seq_pairs[num_seq], // scenes_pairs,             // final QuadCLT [][]           scenes_pairs,
+								mb_max_gain,              // double                       mb_max_gain,
+								batch_mode,               // final boolean                batch_mode,
+								debugLevel);              // final int                    debugLevel)
+						if (sfm_dsn != null) {
+							combo_dsn_final = sfm_dsn;
+							done_sfm = true;
+						}
 					}
 				} 
 				if (!done_sfm) { // first pass or sfm failed
 					// should skip scenes w/o orientation 06/29/2022
 					combo_dsn_final = intersceneExport( // result indexed by COMBO_DSN_TITLES, COMBO_DSN_INDX_***
 							clt_parameters,      // CLTParameters        clt_parameters,
-							ers_reference,       // ErsCorrection        ers_reference,
+//							ers_reference,       // ErsCorrection        ers_reference,
+							ref_index,           // int                  ref_index,							
 							quadCLTs,            // QuadCLT []           scenes,
 							colorProcParameters, // ColorProcParameters  colorProcParameters,
 							debugLevel);         // int                  debug_level
 					if (clt_parameters.imp.sky_recalc) { // force blue sky recalculation even if it exists
-						reuseRefDSI(
+						reuseRefDSI( // processes blue sky
 								clt_parameters,      // CLTParameters  clt_parameters,
 								colorProcParameters, // ColorProcParameters colorProcParameters,
 								quadCLT_main,        // QuadCLT        quadCLT_main, // tiles should be set
@@ -4927,9 +5067,6 @@ public class OpticalFlow {
 								debugLevel);         // int            debugLevel)
 					}
 				}
-				
-				
-				
 				quadCLTs[ref_index].inc_accum();
 				// save with updated num_accum
 				quadCLTs[ref_index].saveInterProperties( // save properties for interscene processing (extrinsics, ers, ...)
@@ -4937,9 +5074,9 @@ public class OpticalFlow {
 						debugLevel+1);
 			}
             double [] reduced_strength = new double[1];
- //   		boolean use_lma_dsi =      clt_parameters.imp.use_lma_dsi;
 			if (quadCLTs[ref_index].getNumOrient() < min_num_orient) {
 				boolean [] reliable_ref = null;
+				boolean show_reliable_ref = false;
 				if (min_ref_str > 0.0) {
 					reliable_ref = quadCLTs[ref_index].getReliableTiles( // will be null if does not exist.
 							use_combo_relaible, // boolean use_combo,
@@ -4947,10 +5084,21 @@ public class OpticalFlow {
 							min_ref_frac,       // double min_ref_frac,
 							true,               // boolean needs_lma);
 							reduced_strength);  // if not null will return >0 if had to reduce strength (no change if did not reduce)
+					
+					if (show_reliable_ref) {
+						double [] dbg_img = new double [reliable_ref.length];
+						for (int i = 0; i < dbg_img.length; i++) {
+							dbg_img[i] = reliable_ref[i]?1:0;
+						}
+						ShowDoubleFloatArrays.showArrays(
+								dbg_img,
+								quadCLTs[ref_index].getTileProcessor().getTilesX(),
+								quadCLTs[ref_index].getTileProcessor().getTilesY(),
+								"reliable_ref-"+quadCLTs[ref_index].getImageName());
+					}
+					
+					
 				}
-//				if (reduced_strength[0] > 0) {
-//					use_lma_dsi = false; // too few points
-//				}
 				// on last pass use final max MB correction same as for render (mb_max_gain - typical =5.0),
 				// for earlier - mb_max_gain_inter (which may be smaller - typical = 2.0)
 				double mb_max_gain = clt_parameters.imp.mb_max_gain;
@@ -4962,6 +5110,8 @@ public class OpticalFlow {
 						mb_max_gain,    // double         mb_max_gain,
 						reliable_ref,   // boolean []     reliable_ref, // null or bitmask of reliable reference tiles
 						quadCLTs,       // QuadCLT []     quadCLTs,
+						ref_index,      // int            ref_index,
+						new int [] {earliest_scene, last_index},  //int     []     range,
 						!batch_mode,    // boolean        test_motion_blur,
 						debugLevel) ;   // int            debugLevel)
 				// should update earliest_scene
@@ -5050,6 +5200,7 @@ public class OpticalFlow {
 					set_channels[ref_index].set_name, // String set_name
 					quadCLTs,                         // QuadCLT[]           quadCLTs,
 					earliest_scene,                   // int                 earliest_scene,
+					last_index,                       // int                 last_scene,
 					ref_index,                        // int                 ref_index,
 					combo_dsn_final_filtered,         // double [][]         combo_dsn_final_filtered,
 					threadsMax,                       // int                 threadsMax,  // maximal number of threads to launch
@@ -5105,11 +5256,11 @@ public class OpticalFlow {
 			test_ers0 = quadCLTs.length -1; // make it always == reference !
 			testERS(
 					clt_parameters, // CLTParameters clt_parameters,
-					test_ers0, // int           indx0, // reference scene in a pair
-					test_ers1, // int           indx1, // other scene in a pair
-					//		    		double []     ref_disparity,			
-					quadCLTs, // QuadCLT []    quadCLTs,
-					debugLevel); // int           debugLevel)
+					test_ers0,      // int           indx0, // reference scene in a pair
+					test_ers1,      // int           indx1, // other scene in a pair
+					quadCLTs,       // QuadCLT []    quadCLTs,
+					ref_index,      // int           ref_index,
+					debugLevel);    // int           debugLevel)
 			System.out.println("buildSeries(): ABORTED after test_ers"); //
 			return quadCLTs[ref_index].getX3dTopDirectory();
 		}
@@ -5128,12 +5279,7 @@ public class OpticalFlow {
 	        	disparity_raw = new double [tilesX * tilesY];
 	        	Arrays.fill(disparity_raw,clt_parameters.disparity);
 	        	if (combo_dsn_final == null) {
-///	        		combo_dsn_final = quadCLTs[ref_index].readDoubleArrayFromModelDirectory(
-///	        				"-INTER-INTRA-LMA", // String      suffix,
-///	        				0, // int         num_slices, // (0 - all)
-///	        				null); // int []      wh);
 					combo_dsn_final =quadCLTs[ref_index].restoreComboDSI(true); // also sets quadCLTs[ref_index].dsi and blue sky
-	        		
 	        	}
 	        	double [][] dls = { 
 	        			combo_dsn_final[COMBO_DSN_INDX_DISP],
@@ -5181,8 +5327,9 @@ public class OpticalFlow {
 	        	// minimum fill peripheral areas with Laplassian?
 	        	double [][] dxyzatr_dt = new double [quadCLTs.length][];
 	        	int [][] min_max_vel = getERSStats(
-	        			clt_parameters,   //  CLTParameters clt_parameters,
+	        			clt_parameters,   // CLTParameters clt_parameters,
 	        			quadCLTs,         // QuadCLT []    quadCLTs,
+	        			ref_index,        // int           ref_index,
 	        			dxyzatr_dt);      // double [][] dxyzatr_dt);
 	        	double [] vel_avg = new double [min_max_vel.length];
 	        	int navg = 0;
@@ -5302,6 +5449,7 @@ public class OpticalFlow {
 	        							scenes_suffix,      // String        suffix,
 	        							ds_vantage[0],      // selected_disparity, // double []     ref_disparity,			
 	        							quadCLTs,           // QuadCLT []    quadCLTs,
+	        					        ref_index,          // int           ref_index,
 	        							threadsMax,         // int           threadsMax,
 	        							debugLevel);        // int           debugLevel);
 	        					if (save_mapped_mono_color[col_mode]) {	        	
@@ -5627,9 +5775,10 @@ public class OpticalFlow {
 					clt_parameters,      // CLTParameters                            clt_parameters,
 					colorProcParameters, // ColorProcParameters                      colorProcParameters,
 					rgbParameters,       // EyesisCorrectionParameters.RGBParameters rgbParameters,
-					quadCLTs[quadCLTs.length-1], // quadCLT_main,        // final QuadCLT                            parameter_scene, // to use for rendering parameters in multi-series sequences
+					quadCLTs[ref_index], // quadCLTs.length-1], // quadCLT_main,        // final QuadCLT                            parameter_scene, // to use for rendering parameters in multi-series sequences
 		            // if null - use reference scene 
 					quadCLTs,            // QuadCLT []                               scenes,
+					ref_index,           // final int                                ref_index,
 					combo_dsn_final,     //double [][]                              combo_dsn_final, // null OK, will read file
 					updateStatus,        // final boolean                            updateStatus,
 					debugLevel); //  + 1);         // final int                                debugLevel)
@@ -6301,6 +6450,7 @@ public class OpticalFlow {
     		int           indx0, // reference scene in a pair
     		int           indx1, // other scene in a pair
     		QuadCLT []    quadCLTs,
+    		int           ref_index,
     		int           debugLevel) {
     	// First create a pair of images, similar to renderSceneSequence()
 		boolean show_color =  clt_parameters.imp.show_mapped_color;
@@ -6308,7 +6458,6 @@ public class OpticalFlow {
 	    boolean use_combo_dsi =        clt_parameters.imp.use_combo_dsi;
 	    boolean use_lma_dsi =          clt_parameters.imp.use_lma_dsi;
 		
-    	int ref_index = quadCLTs.length-1;
     	indx0 = ref_index; // disregard initial setting, set to reference
 		int tilesX =  quadCLTs[ref_index].getTileProcessor().getTilesX();
         int tilesY =  quadCLTs[ref_index].getTileProcessor().getTilesY();
@@ -6591,8 +6740,8 @@ public class OpticalFlow {
     public static int [][] getERSStats(
     		CLTParameters clt_parameters,
     		QuadCLT []    quadCLTs,
-    		double [][] dxyzatr_dt) {
-        int ref_index = quadCLTs.length -1;
+    		int           ref_index,
+    		double [][]   dxyzatr_dt) {
         int earlies_scene = quadCLTs[ref_index].getEarliestScene(quadCLTs);
 		ErsCorrection ers_reference = quadCLTs[ref_index].getErsCorrection();
         int dbg_scene = -95;
@@ -6655,6 +6804,7 @@ public class OpticalFlow {
     		String         suffix_in,
     		double []      ref_disparity,			
     		QuadCLT []     quadCLTs,
+    		int            ref_index,
     		int            threadsMax,
     		int            debugLevel) {
     	boolean corr_raw_ers = true;
@@ -6679,7 +6829,6 @@ public class OpticalFlow {
     	if ((mode3d <0) && (corr_raw_ers)) {
     		suffix+="-RAW_ERS";
     	}
-        int ref_index = quadCLTs.length -1;
     	int num_sens = quadCLTs[ref_index].getNumSensors();
 		ErsCorrection ers_reference = quadCLTs[ref_index].getErsCorrection();
         int num_used_sens = 0;
@@ -8160,8 +8309,9 @@ public class OpticalFlow {
 	
 	public static double [][] intersceneExport(
 			CLTParameters        clt_parameters,
-			ErsCorrection        ers_reference,
-			QuadCLT []           scenes,
+//			ErsCorrection        ers_reference,
+			int                  ref_index,
+			QuadCLT []           scenes_in,
 			ColorProcParameters  colorProcParameters,
 			int                  debug_level
 			)
@@ -8186,13 +8336,24 @@ public class OpticalFlow {
 		case 1:   fsplit_refines = 1; break;
 		case 2:   fsplit_refines = 3; break;
 		}
-		
+		QuadCLT []           scenes = scenes_in.clone();
+		ErsCorrection        ers_reference = scenes[ref_index].getErsCorrection();
+		int      indx_ref = scenes.length - 1; // Always added to the end even if out-of order
+		if (ref_index != indx_ref) {
+			System.arraycopy(
+					scenes_in,
+					ref_index+1,
+					scenes,
+					ref_index,
+					indx_ref-ref_index);
+			scenes[indx_ref] = scenes_in[ref_index];
+		}
 		final double scale_far_split_strength = 1.5;
 		// empirical correction for both lma and non-lma step
 	  	double corr_nonlma = 1.0; // 1.23;
 	  	double corr_lma =    1.0; // 1.23;
 		// reference scene is always added to the end, even is out of timestamp order
-		int      indx_ref = scenes.length - 1; // Always added to the end even if out-of order
+		
 		
 		QuadCLT  ref_scene = scenes[indx_ref]; // ordered by increasing timestamps
 		boolean generate_outlines = false; // true; // TODO: move to configs
@@ -8333,7 +8494,9 @@ public class OpticalFlow {
 						combo_dsn_titles_full.length:combo_dsn_titles.length][combo_dsn[0].length];
 		combo_dsn_final[COMBO_DSN_INDX_DISP]= combo_dsn[COMBO_DSN_INDX_DISP].clone();
 		for (int i = 1; i < combo_dsn_final.length; i++) {
-			Arrays.fill(combo_dsn_final[i], Double.NaN);
+			if (i != COMBO_DSN_INDX_SFM_GAIN) {
+				Arrays.fill(combo_dsn_final[i], Double.NaN);
+			}
 		}
 		// Save pair selection and minimize them for scanning, then restore;
 		int num_sensors =scenes[indx_ref].getNumSensors();

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

@@ -186,7 +186,106 @@ public class QuadCLTCPU {
     public Did_gps_pos did_gps1_ubx_pos = null;
     public String      ims_last_path =    null;
     
+    // latest (or any? first also?) scene in a sequence has a timestamp of the reference scene  
+    public String      timestamp_reference = null;
+    // only reference scene has a pair of first/last scene in a sequence
+    public String      timestamp_first =     null;
+    public String      timestamp_last =      null;
     
+    public String getReferenceTimestamp() {
+    	return timestamp_reference;
+    }
+    
+    /**
+     * If this scene has a reference to one of the scenes array, return its index in the array
+     * otherwise return -1 if there is no pointer to reference scene and -2 if pointer is not
+     * within the provided array. OK to call with scenes==null, results can be -1 (no link) and -2
+     * @param scenes array of scenes (QuadCLTCPU instances)
+     * @return index of the reference scene in the provided scenes array or negative on error
+     */
+    public int getReferenceIndex(QuadCLTCPU [] scenes) {
+    	if (timestamp_reference != null) {
+    		if (scenes != null) {
+    			for (int i = 0; i < scenes.length; i++) {
+    				if (timestamp_reference.equals(scenes[i].getImageName())) {
+    					return i;
+    				}
+    			}
+    		}
+    		return -2;
+    	}
+    	return -1;
+    }
+    
+    /**
+     * If this scene is a reference scene, then return array of {first_index, last_index}.
+     * If this scene is not a reference or the array does not contain both the first and
+     * the last scene - return null.
+     * @param scenes array of scenes (QuadCLTCPU instances)
+     * @return a pair of {first_index, last_index}
+     */
+    public int [] getFirstLastIndex(QuadCLTCPU [] scenes) { // should be ordered accending ts
+		int [] fl = new int [] {-1,-1};
+    	if ((timestamp_first != null) && (timestamp_last != null)) {
+    		int i = 0;
+    		for (; i < scenes.length; i++) {
+    			if (timestamp_first.equals(scenes[i].getImageName())) {
+    				fl[0] = i;
+    				break;
+    			}
+    		}
+    		for (; i < scenes.length; i++) {
+    			if (timestamp_last.equals(scenes[i].getImageName())) {
+    				fl[1] = i;
+    				break;
+    			}
+    		}
+    	}
+    	if ((fl[0] >= 0) && (fl[1] >= 0)) {
+    		return fl;
+    	} else {
+    		return null;
+    	}
+    }
+
+    /**
+     * Set reference scene pointer by a reference scene timestamp.
+     * @param ts reference scene pointer as String
+     */
+    public void setRefPointer(String ts) {
+    	timestamp_reference = ts;
+    }
+    
+    /**
+     * Set reference scene pointer by a reference scene instance.
+     * @param scene reference scene instance
+     */
+    public void setRefPointer(QuadCLTCPU scene) {
+    	setRefPointer(scene.getImageName());
+    }
+    
+    /**
+     * Set first, last scene pointers by timestamps.
+     * @param ts0 first scene pointer as String
+     * @param ts1 last scene pointer as String
+     */
+    public void setFirstLastPointers(String ts0,String ts1) {
+    	timestamp_first = ts0;
+    	timestamp_last = ts1;
+    }
+
+    /**
+     * Set first, last scene pointers by instances.
+     * @param scene0 first scene instance
+     * @param scene1 last scene instance
+     */
+    public void setFirstLastPointers(QuadCLTCPU scene0,QuadCLTCPU scene1) {
+    	setFirstLastPointers(
+    			scene0.getImageName(),
+    			scene1.getImageName());
+    }
+    
+    //QuadCLTCPU [] scenes
     // TODO: Use dsi[] instead
 //    @Deprecated
 //    public boolean[] blue_sky =                            null;
@@ -1666,21 +1765,22 @@ public class QuadCLTCPU {
 				}
 			}
 		}
-		int num_sfm_gain = 0;
-		for (int i = 0; i < reliable.length; i++) {
-			if (sfm_gain[i] >0) {
-				num_sfm_gain++;
-			}
-		}
-		if (num_sfm_gain > min_reliable) {
+		if (sfm_gain != null) {
+			int num_sfm_gain = 0;
 			for (int i = 0; i < reliable.length; i++) {
-				if (sfm_gain[i] <= 0){
-					reliable[i] = false;
-					strength[i] = 0.0;
+				if (sfm_gain[i] >0) {
+					num_sfm_gain++;
+				}
+			}
+			if (num_sfm_gain > min_reliable) {
+				for (int i = 0; i < reliable.length; i++) {
+					if (sfm_gain[i] <= 0){
+						reliable[i] = false;
+						strength[i] = 0.0;
+					}
 				}
 			}
 		}
-		
 		int num_reliable = 0;
 		for (boolean b: reliable) if (b) num_reliable++;
 		if (num_reliable < min_reliable) { // not enough, select best tiles, ignoring LMA
@@ -1859,13 +1959,10 @@ public class QuadCLTCPU {
 	public Properties restoreInterProperties( // restore properties for interscene processing (extrinsics, ers, ...)
 			String path,             // full name with extension or null to use x3d directory
 			boolean all_properties,
-			
 			int debugLevel)
 	{
 		if (path == null) {
-//			path = image_name + ((properties == null) ? "-INTERFRAME":"")+".corr-xml";
 			path = image_name + (all_properties? "": "-INTERFRAME")+".corr-xml";
-
 		}
 		if (!path.contains(Prefs.getFileSeparator())) {
 			String x3d_path = getX3dDirectory();
@@ -1874,11 +1971,11 @@ public class QuadCLTCPU {
 		properties = loadProperties(
 				path, // String path,
 				properties, // Properties properties)
-				false); // boolean silent)
-
-//		if (properties == null) {
-//			properties = new Properties();
-//		}
+				true); // false); // boolean silent)
+		if (properties == null ) {
+			System.out.println("Failed to restore interframe properties from :"+path);
+			return null;
+		}
 		String prefix = is_aux?PREFIX_AUX:PREFIX; 
 		getProperties(prefix); // will set Geometry correction non-null
 		if (!(geometryCorrection instanceof ErsCorrection)) { // should only be for the new GeometryCorrection created in getProperties
@@ -3252,6 +3349,16 @@ public class QuadCLTCPU {
 		if (this.photometric_scene != null) {
 			properties.setProperty(prefix+"photometric_scene", this.photometric_scene);
 		}
+		
+		if (this.timestamp_reference != null) {
+			properties.setProperty(prefix+"timestamp_reference", this.timestamp_reference);
+		}
+		if (this.timestamp_first != null) {
+			properties.setProperty(prefix+"timestamp_first", this.timestamp_first);
+		}
+		if (this.timestamp_last != null) {
+			properties.setProperty(prefix+"timestamp_last", this.timestamp_last);
+		}
 		return properties;
 	}
 
@@ -3406,6 +3513,16 @@ public class QuadCLTCPU {
 		if (properties.getProperty(prefix+"photometric_scene")!=null) {
 			this.photometric_scene= (String) properties.getProperty(prefix+"photometric_scene");
 		}
+		if (properties.getProperty(prefix+"timestamp_reference")!=null) {
+			this.timestamp_reference= (String) properties.getProperty(prefix+"timestamp_reference");
+		}
+		if (properties.getProperty(prefix+"timestamp_first")!=null) {
+			this.timestamp_first= (String) properties.getProperty(prefix+"timestamp_first");
+		}
+		if (properties.getProperty(prefix+"timestamp_last")!=null) {
+			this.timestamp_last= (String) properties.getProperty(prefix+"timestamp_last");
+		}
+
 	}
 
 	public void setKernelImageFile(ImagePlus img_kernels){ // not used in lwir
@@ -9526,9 +9643,14 @@ public class QuadCLTCPU {
 
 // need to replace for low-res?
 				  preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						  imp_srcs, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						  saturation_imp, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						  clt_parameters,
+						  // adding these parameters for more flexibility in accuracy/speed
+						  clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						  clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						  clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						  clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						  false,                              // boolean no_bg_generate,			  
+						  false,                              // boolean no_lma,			  
 						  debayerParameters,
 						  colorProcParameters,
 						  rgbParameters,
@@ -9826,6 +9948,13 @@ public class QuadCLTCPU {
 
 	  public boolean preExpandCLTQuad3d( // USED in lwir
 			  CLTParameters                                    clt_parameters,
+			  // adding these parameters for more flexibility in accuracy/speed
+			  int                                              gr_max_clust_radius,			  
+			  double                                           disp_scan_start,
+			  double                                           disp_scan_step,
+			  int                                              disp_scan_count,
+			  boolean                                          no_bg_generate,
+			  boolean                                          no_lma,
 			  EyesisCorrectionParameters.DebayerParameters     debayerParameters,
 			  ColorProcParameters                              colorProcParameters,
 			  EyesisCorrectionParameters.RGBParameters         rgbParameters,
@@ -9847,24 +9976,28 @@ public class QuadCLTCPU {
 		  final double     fom_cdiff =         clt_parameters.getFomCDiff(isLwir()); // 0.02;
 		  final double     fom_inf_bonus  =    clt_parameters.fom_inf_bonus;  // 0.2; // add this to infinity FOM (if it is closer than fom_inf_range)
 		  final double     fom_inf_range  =    clt_parameters.fom_inf_range;  // 0.5;
-
+		  
+//		  final boolean no_lma = true;
 		  //max_expand
 //		  String name = (String) imp_quad[0].getProperty("name");
 		  String name = getImageName();
 		  // should create data for the macro! (diff, rgb) make sure .texture_tiles is measured correctly
-		  CLTPass3d bgnd_data_lma = CLTBackgroundMeas( // measure background - both CPU and GPU (remove textures from GPU)
-				  clt_parameters,
-				  true,              //  final boolean     run_lma, //			  
-				  max_mismatch,      // max_chn_diff, //final double      max_chn_diff, // filter correlation results by maximum difference between channels
-				  mismatch_override, // final double      mismatch_override, // keep tile with large mismatch if there is LMA with really strong correlation
-				  threadsMax,        // maximal number of threads to launch
-				  updateStatus,
-				  debugLevel);
-		  if (clt_parameters.img_dtt.lmamask_dbg) {
-			  System.out.println("Remove me - QCC8257");
-			  return false;
+		  CLTPass3d bgnd_data_lma = null;
+		  if (!no_lma) {
+			  bgnd_data_lma = CLTBackgroundMeas( // measure background - both CPU and GPU (remove textures from GPU)
+					  clt_parameters,
+					  true,              //  final boolean     run_lma, //			  
+					  max_mismatch,      // max_chn_diff, //final double      max_chn_diff, // filter correlation results by maximum difference between channels
+					  mismatch_override, // final double      mismatch_override, // keep tile with large mismatch if there is LMA with really strong correlation
+					  threadsMax,        // maximal number of threads to launch
+					  updateStatus,
+					  debugLevel);
+			  if (clt_parameters.img_dtt.lmamask_dbg) {
+				  System.out.println("Remove me - QCC8257");
+				  return false;
+			  }
+			  tp.clt_3d_passes.add(bgnd_data_lma);
 		  }
-		  tp.clt_3d_passes.add(bgnd_data_lma);
 		  CLTPass3d bgnd_data_nolma = CLTBackgroundMeas( // measure background - both CPU and GPU (remove textures from GPU)
 				  clt_parameters,
 				  false,              //  final boolean     run_lma, //			  
@@ -9874,102 +10007,95 @@ public class QuadCLTCPU {
 				  updateStatus,
 				  debugLevel);
 		  tp.clt_3d_passes.add(bgnd_data_nolma);
-		  
-		  //    	  if (show_init_refine)
-//		  if ((debugLevel > -2) && clt_parameters.show_first_bg) {
-		  if ((debugLevel > -2) && clt_parameters.show_first_bg) {
-			  tp.showScan(
-					  tp.clt_3d_passes.get(0), // CLTPass3d   scan,
-					  "bgnd_data_lma-"+tp.clt_3d_passes.size());
-			  tp.showScan(
-					  tp.clt_3d_passes.get(1), // CLTPass3d   scan,
-					  "bgnd_data_nolma-"+tp.clt_3d_passes.size());
+		  if (no_lma) {
+			  bgnd_data_lma = bgnd_data_nolma; // use for bg generation
 		  }
+		  if (!no_bg_generate) {
+			  //    	  if (show_init_refine)
+			  //		  if ((debugLevel > -2) && clt_parameters.show_first_bg) {
+			  if ((debugLevel > -2) && clt_parameters.show_first_bg) {
+				  tp.showScan(
+						  tp.clt_3d_passes.get(0), // CLTPass3d   scan,
+						  "bgnd_data_lma-"+tp.clt_3d_passes.size());
+				  tp.showScan(
+						  tp.clt_3d_passes.get(1), // CLTPass3d   scan,
+						  "bgnd_data_nolma-"+tp.clt_3d_passes.size());
+			  }
 
-		  //TODO: Move away from here?
-		  boolean no_image_save = false; // Restore? true;
+			  //TODO: Move away from here?
+			  boolean no_image_save = false; // Restore? true;
 
-		  boolean [] bgmask = getBackgroundImageMasks(
-    			  clt_parameters,
-    			  bgnd_data_lma,   // CLTPass3d               bgnd_data,
-    			  name,        // .getTitle(), //String name=(String) imp_src.getProperty("name");
-    			  threadsMax,  // maximal number of threads to launch
-    			  updateStatus,
-    			  debugLevel);
-		  
-		  bgnd_data_nolma.setSelected(bgnd_data_lma.getSelected());
-		  
-		  if ((debugLevel > -2)  && (bgmask != null)) {
-			  double [][] dbg_img = new double[1][tp.getTilesY() * tp.getTilesX()];
-			  for (int i = 0; i<dbg_img[0].length;i++){
-				  dbg_img[0][i] =  bgmask[i]?1:0;
+			  boolean [] bgmask = getBackgroundImageMasks(
+					  clt_parameters,
+					  bgnd_data_lma,   // CLTPass3d               bgnd_data,
+					  name,        // .getTitle(), //String name=(String) imp_src.getProperty("name");
+					  threadsMax,  // maximal number of threads to launch
+					  updateStatus,
+					  debugLevel);
+
+			  bgnd_data_nolma.setSelected(bgnd_data_lma.getSelected());
+
+			  if ((debugLevel > -2)  && (bgmask != null)) {
+				  double [][] dbg_img = new double[1][tp.getTilesY() * tp.getTilesX()];
+				  for (int i = 0; i<dbg_img[0].length;i++){
+					  dbg_img[0][i] =  bgmask[i]?1:0;
+				  }
+				  ShowDoubleFloatArrays.showArrays(
+						  dbg_img,
+						  tp.getTilesX(),
+						  tp.getTilesY(),
+						  true, "bgmask");
 			  }
-			  ShowDoubleFloatArrays.showArrays(
-					dbg_img,
-					tp.getTilesX(),
-					tp.getTilesY(),
-					true, "bgmask");
-		  }
-		  
-    	  ImagePlus imp_bgnd_int = getBackgroundImage( // null pointer
-    			  bgmask, // boolean []                                bgnd_tiles, 
-    			  bgnd_data_lma,   // CLTPass3d               bgnd_data,
-    			  clt_parameters,
-    			  colorProcParameters,
-    			  rgbParameters,
-    			  name,               // .getTitle(), //String name=(String) imp_src.getProperty("name");
-    			  threadsMax,  // maximal number of threads to launch
-    			  updateStatus,
-    			  debugLevel);
-		  
-    	  if (debugLevel > -2) {
-    		  imp_bgnd_int.show(); /// OK 
-    	  }
-		  
-//		  if (debugLevel > -100) {
-//			  return null;
-//		  }
-		  // resize for backdrop here! check imp_bgnd_int here !
 
-    	  ImagePlus imp_bgnd = finalizeBackgroundImage( // USED in lwir - output all 0?
-    			  imp_bgnd_int,   //  ImagePlus imp_texture_bgnd, null pointer
-    			  no_image_save,  // boolean    no_image_save,
-    			  clt_parameters, // CLTParameters           clt_parameters,
-    			  name,           // String     name,
-    			  debugLevel);    // int        debugLevel)
+			  ImagePlus imp_bgnd_int = getBackgroundImage( // null pointer
+					  bgmask, // boolean []                                bgnd_tiles, 
+					  bgnd_data_lma,   // CLTPass3d               bgnd_data,
+					  clt_parameters,
+					  colorProcParameters,
+					  rgbParameters,
+					  name,               // .getTitle(), //String name=(String) imp_src.getProperty("name");
+					  threadsMax,  // maximal number of threads to launch
+					  updateStatus,
+					  debugLevel);
 
-////    	  imp_bgnd.show(); /// all black -> change clt_parameters.alpha0, clt_parameters.alpha1 
-    	  
-    	  bgnd_data_lma.texture = (imp_bgnd == null)? null: ( imp_bgnd.getTitle()+ (clt_parameters.black_back? ".jpeg" : ".png"));
+			  if (debugLevel > -2) {
+				  imp_bgnd_int.show(); /// OK 
+			  }
 
-    	  // create x3d file
-    	  X3dOutput x3dOutput = new X3dOutput(
-    			  clt_parameters,
-    			  correctionsParameters,
-    			  geometryCorrection,
-    			  tp.clt_3d_passes);
+			  //		  if (debugLevel > -100) {
+			  //			  return null;
+			  //		  }
+			  // resize for backdrop here! check imp_bgnd_int here !
 
-		  x3dOutput.generateBackground(clt_parameters.infinityDistance <= 0.0); // needs just first (background) scan
+			  ImagePlus imp_bgnd = finalizeBackgroundImage( // USED in lwir - output all 0?
+					  imp_bgnd_int,   //  ImagePlus imp_texture_bgnd, null pointer
+					  no_image_save,  // boolean    no_image_save,
+					  clt_parameters, // CLTParameters           clt_parameters,
+					  name,           // String     name,
+					  debugLevel);    // int        debugLevel)
 
+			  bgnd_data_lma.texture = (imp_bgnd == null)? null: ( imp_bgnd.getTitle()+ (clt_parameters.black_back? ".jpeg" : ".png"));
 
+		  }
+
+		  // create x3d file
+		  X3dOutput x3dOutput = new X3dOutput(
+				  clt_parameters,
+				  correctionsParameters,
+				  geometryCorrection,
+				  tp.clt_3d_passes);
+
+		  x3dOutput.generateBackground(clt_parameters.infinityDistance <= 0.0); // needs just first (background) scan
 
     	  // refine first measurement
     	  int bg_pass = tp.clt_3d_passes.size() - 1; // 0
     	  int refine_pass = tp.clt_3d_passes.size(); // 1
-
-    	  //		  final boolean show_init_refine = true;
-    	  //		  final boolean show_expand =      true;
-
-    	  //    	  if (show_init_refine)
     	  if (debugLevel > -1) {
     		  tp.showScan(
     				  tp.clt_3d_passes.get(bg_pass), // CLTPass3d   scan,
     				  "after_bg-"+tp.clt_3d_passes.size());
     	  }
 
-//    	  final double     weight_var = 1.0; // 1.0;   // weight of variance data (old, detects thin wires?)
- //   	  final double     weight_Y =   1.0;     // weight of average intensity
- //   	  final double     weight_RBmG = 5.0;  // weight of average color difference (0.5*(R+B)-G), shoukld be ~5*weight_Y
 
     	  // TODO: Make double pass - with only 	weight_var (thin wires) and weight_Y, weight_RBmG - larger objects
     	  // just use two instances of MacroCorrelation, run one after another (move code to MacroCorrelation class)
@@ -9977,9 +10103,9 @@ public class QuadCLTCPU {
     	  ArrayList <CLTPass3d> new_meas = null; // filled either from macro correlation, or just from plain disparity scan
     	  if (no_macro) {
     		  new_meas = prepareDisparityScan(
-    				  clt_parameters.disp_scan_start, // double scan_start,
-    				  clt_parameters.disp_scan_step, // 	double scan_step,
-    				  clt_parameters.disp_scan_count); // 	int    scan_count)
+    				  disp_scan_start,  // double scan_start,
+    				  disp_scan_step,   // 	double scan_step,
+    				  disp_scan_count); // 	int    scan_count)
     	  } else {
     		  MacroCorrelation mc = new MacroCorrelation(
     				  tp,
@@ -10292,8 +10418,6 @@ public class QuadCLTCPU {
     			  break;
     		  }
     	  }
-
-    	  
     	  
 /// Refining after all added
 		  if (debugLevel > -3){
@@ -10550,11 +10674,10 @@ public class QuadCLTCPU {
 		  
 		  
 //		  int max_clust_radius =  4; // 7x7
-		  CLTPass3d [] combo_multi = new CLTPass3d[clt_parameters.gr_max_clust_radius+1];
+		  CLTPass3d [] combo_multi = new CLTPass3d[gr_max_clust_radius+1];
 		  combo_multi[0] = combo_pass;
 		  int max_expand_radius = 0; // max_clust_radius;
-		  for (int clust_radius = 2; clust_radius <= clt_parameters.gr_max_clust_radius; clust_radius++) {
-
+		  for (int clust_radius = 2; clust_radius <= gr_max_clust_radius; clust_radius++) {
 	    	  if (clt_parameters.gr_reduce_multi) {
 		    	  clt_parameters.img_dtt.setMcorr(getNumSensors(), 0 ); // remove all
 		    	  clt_parameters.img_dtt.setMcorrNeib(getNumSensors(),true);
@@ -10562,7 +10685,6 @@ public class QuadCLTCPU {
 		    	  clt_parameters.img_dtt.setMcorrDia (getNumSensors(),true); // remove even more?
 	    	  }
 			  
-			  
 			  // using combo_pass (latest)
 			  int num_added = 0;
 			  if (clust_radius < max_expand_radius) { // disabled currently by int max_expand_radius = 0;
@@ -15337,9 +15459,14 @@ public class QuadCLTCPU {
 			  if (correctionsParameters.clt_batch_extrinsic) {
 				  if (tp != null) tp.resetCLTPasses();
 				  boolean ok = preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						  imp_srcs, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						  saturation_imp, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						  clt_parameters,
+						  // adding these parameters for more flexibility in accuracy/speed
+						  clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						  clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						  clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						  clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						  false,                              // boolean no_bg_generate,			  
+						  false,                              // boolean no_lma,			  
 						  debayerParameters,
 						  colorProcParameters,
 						  rgbParameters,
@@ -15361,9 +15488,14 @@ public class QuadCLTCPU {
 			  if (correctionsParameters.clt_batch_poly) {
 				  if (tp != null) tp.resetCLTPasses();
 				  boolean ok = preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						  imp_srcs, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						  saturation_imp, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						  clt_parameters,
+						  // adding these parameters for more flexibility in accuracy/speed
+						  clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						  clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						  clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						  clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						  false,                              // boolean no_bg_generate,			  
+						  false,                              // boolean no_lma,			  
 						  debayerParameters,
 						  colorProcParameters,
 						  rgbParameters,
@@ -15402,9 +15534,14 @@ public class QuadCLTCPU {
 			  if (correctionsParameters.clt_batch_explore) {
 				  if (tp != null) tp.resetCLTPasses();
 				  boolean ok = preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						  imp_srcs, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						  saturation_imp, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						  clt_parameters,
+						  // adding these parameters for more flexibility in accuracy/speed
+						  clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						  clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						  clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						  clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						  false,                              // boolean no_bg_generate,			  
+						  false,                              // boolean no_lma,			  
 						  debayerParameters,
 						  colorProcParameters,
 						  rgbParameters,

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

@@ -2447,6 +2447,7 @@ public class TexturedModel {
 			final QuadCLT                            parameter_scene, // to use for rendering parameters in multi-series sequences
 			// if null - use reference scene 
 			QuadCLT []                               scenes,
+			final int                                ref_index,
 			double [][]                              combo_dsn_final, // null OK, will read file
 			final boolean                            updateStatus,
 			final int                                debugLevel)
@@ -2512,7 +2513,7 @@ public class TexturedModel {
 		final boolean tex_um_fixed =        clt_parameters.tex_um_fixed;  // imp.mono_fixed; //  true; // normalize to fixed range when converting to 8 bits 
 		final double  tex_um_range =        clt_parameters.tex_um_range;  // imp.mono_range; // 500.0;  // monochrome full-scale range (+/- half)
 
-		final int           ref_index =  scenes.length - 1;
+//		final int           ref_index =  scenes.length - 1;
 		final QuadCLT       ref_scene =  scenes[ref_index];
 		final TileProcessor tp =         ref_scene.getTileProcessor();
 

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

@@ -7597,9 +7597,14 @@ if (debugLevel > -100) return true; // temporarily !
 				}
 
 				quadCLT_main.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						imp_srcs_main, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						saturation_imp_main, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						clt_parameters,
+						// adding these parameters for more flexibility in accuracy/speed
+						clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						false,                              // boolean no_bg_generate,			  
+						false,                              // boolean no_lma,			  
 						debayerParameters,
 						colorProcParameters,
 						rgbParameters,
@@ -7636,9 +7641,14 @@ if (debugLevel > -100) return true; // temporarily !
 							" (w/o rig), pass "+(num_adjust_aux+1)+" of "+quadCLT_main.correctionsParameters.rig_batch_adjust_aux);
 				}
 				quadCLT_aux.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						imp_srcs_aux, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						saturation_imp_aux, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						clt_parameters,
+						// adding these parameters for more flexibility in accuracy/speed
+						clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						false,                              // boolean no_bg_generate,			  
+						false,                              // boolean no_lma,			  
 						debayerParameters,
 						colorProcParameters_aux,
 						rgbParameters,
@@ -7714,9 +7724,14 @@ if (debugLevel > -100) return true; // temporarily !
 				{
 					if (updateStatus) IJ.showStatus("Building basic DSI for the main camera image set "+quadCLT_main.image_name+ " (post rig adjustment)");
 					quadCLT_main.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//							imp_srcs_main, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//							saturation_imp_main, // boolean [][] saturation_imp, // (near) saturated pixels or null
 							clt_parameters,
+							// adding these parameters for more flexibility in accuracy/speed
+							clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+							clt_parameters.disp_scan_start,     // double  disp_scan_start,
+							clt_parameters.disp_scan_step,      // double  disp_scan_step,
+							clt_parameters.disp_scan_count,     // double  disp_scan_count,
+							false,                              // boolean no_bg_generate,			  
+							false,                              // boolean no_lma,			  
 							debayerParameters,
 							colorProcParameters,
 							rgbParameters,
@@ -7874,9 +7889,14 @@ if (debugLevel > -100) return true; // temporarily !
 			if (quadCLT_main.correctionsParameters.clt_batch_dsi_aux) {
 				if (updateStatus) IJ.showStatus("Building basic DSI for the aux camera image set "+quadCLT_main.image_name+" (for DSI export)");
 				quadCLT_aux.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						imp_srcs_aux, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						saturation_imp_aux, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						clt_parameters,
+						// adding these parameters for more flexibility in accuracy/speed
+						clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						false,                              // boolean no_bg_generate,			  
+						false,                              // boolean no_lma,			  
 						debayerParameters,
 						colorProcParameters_aux,
 						rgbParameters,
@@ -7905,9 +7925,14 @@ if (debugLevel > -100) return true; // temporarily !
 			if (quadCLT_main.correctionsParameters.clt_batch_explore) {
 				if (updateStatus) IJ.showStatus("Building basic DSI for the main camera image set "+quadCLT_main.image_name+" (after all adjustments)");
 				quadCLT_main.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//						imp_srcs_main, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//						saturation_imp_main, // boolean [][] saturation_imp, // (near) saturated pixels or null
 						clt_parameters,
+						// adding these parameters for more flexibility in accuracy/speed
+						clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						false,                              // boolean no_bg_generate,			  
+						false,                              // boolean no_lma,			  
 						debayerParameters,
 						colorProcParameters,
 						rgbParameters,
@@ -9335,7 +9360,8 @@ if (debugLevel > -100) return true; // temporarily !
 			// Does not need to specify reference scene - it is always the last even if out of timestamp order
 			OpticalFlow.intersceneExport(
 					clt_parameters,            // CLTParameters       clt_parameters,
-					ers_reference,             // ErsCorrection        ers_reference,
+//					ers_reference,       // ErsCorrection        ers_reference,
+					indx_ref, // ref_index,           // int                  ref_index,							
 					scenes,                    // QuadCLT []           scenes,
 					colorProcParameters,       // ColorProcParameters colorProcParameters,
 					clt_parameters.inp.noise_debug_level // clt_parameters.ofp.debug_level_optical - 1); // 1); // -1); // int debug_level);
@@ -11941,9 +11967,14 @@ if (debugLevel > -100) return true; // temporarily !
 						}
 						//Generates background image in model tree - should be done later, after adjustment (It is overwritten later, so OK)
 						quadCLT_main.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//								imp_srcs_main, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//								saturation_imp_main, // boolean [][] saturation_imp, // (near) saturated pixels or null
 								clt_parameters,
+								// adding these parameters for more flexibility in accuracy/speed
+								clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+								clt_parameters.disp_scan_start,     // double  disp_scan_start,
+								clt_parameters.disp_scan_step,      // double  disp_scan_step,
+								clt_parameters.disp_scan_count,     // double  disp_scan_count,
+								false,                              // boolean no_bg_generate,			  
+								false,                              // boolean no_lma,			  
 								debayerParameters,
 								colorProcParameters,
 								rgbParameters,
@@ -11993,9 +12024,14 @@ if (debugLevel > -100) return true; // temporarily !
 					if (quadCLT_main.correctionsParameters.clt_batch_dsi1){
 						System.out.println("Trying experimental features DSI/ERS");
 						quadCLT_main.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//								imp_srcs_main, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//								saturation_imp_main, // boolean [][] saturation_imp, // (near) saturated pixels or null
 								clt_parameters,
+								// adding these parameters for more flexibility in accuracy/speed
+								clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+								clt_parameters.disp_scan_start,     // double  disp_scan_start,
+								clt_parameters.disp_scan_step,      // double  disp_scan_step,
+								clt_parameters.disp_scan_count,     // double  disp_scan_count,
+								false,                              // boolean no_bg_generate,			  
+								false,                              // boolean no_lma,			  
 								debayerParameters,
 								colorProcParameters,
 								rgbParameters,
@@ -12074,9 +12110,14 @@ if (debugLevel > -100) return true; // temporarily !
 					if (quadCLT_main.correctionsParameters.clt_batch_explore) {
 						if (updateStatus) IJ.showStatus("Building basic DSI for the main camera image set "+quadCLT_main.image_name+" (after all adjustments)");
 						quadCLT_main.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//								imp_srcs_main, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//								saturation_imp_main, // boolean [][] saturation_imp, // (near) saturated pixels or null
 								clt_parameters,
+								// adding these parameters for more flexibility in accuracy/speed
+								clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+								clt_parameters.disp_scan_start,     // double  disp_scan_start,
+								clt_parameters.disp_scan_step,      // double  disp_scan_step,
+								clt_parameters.disp_scan_count,     // double  disp_scan_count,
+								false,                              // boolean no_bg_generate,			  
+								false,                              // boolean no_lma,			  
 								debayerParameters,
 								colorProcParameters,
 								rgbParameters,
@@ -12280,9 +12321,14 @@ if (debugLevel > -100) return true; // temporarily !
 						//			    dbg_path = "/home/elphel/lwir16-proc/proc1/results_cuda/25/extrinsics_bgnd_combo.tif";
 						if (dbg_path == null) {
 							quadCLT_aux.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//									imp_srcs_aux, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//									saturation_imp_aux, // boolean [][] saturation_imp, // (near) saturated pixels or null
 									clt_parameters,
+									// adding these parameters for more flexibility in accuracy/speed
+									clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+									clt_parameters.disp_scan_start,     // double  disp_scan_start,
+									clt_parameters.disp_scan_step,      // double  disp_scan_step,
+									clt_parameters.disp_scan_count,     // double  disp_scan_count,
+									false,                              // boolean no_bg_generate,			  
+									false,                              // boolean no_lma,			  
 									debayerParameters,
 									colorProcParameters_aux,
 									rgbParameters,
@@ -12372,9 +12418,14 @@ if (debugLevel > -100) return true; // temporarily !
 					if (quadCLT_main.correctionsParameters.clt_batch_dsi_aux) {
 						if (updateStatus) IJ.showStatus("Building basic DSI for the aux camera image set "+quadCLT_main.image_name+" (for DSI export)");
 						quadCLT_aux.preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
-//								imp_srcs_aux, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
-//								saturation_imp_aux, // boolean [][] saturation_imp, // (near) saturated pixels or null
 								clt_parameters,
+								// adding these parameters for more flexibility in accuracy/speed
+								clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+								clt_parameters.disp_scan_start,     // double  disp_scan_start,
+								clt_parameters.disp_scan_step,      // double  disp_scan_step,
+								clt_parameters.disp_scan_count,     // double  disp_scan_count,
+								false,                              // boolean no_bg_generate,			  
+								false,                              // boolean no_lma,			  
 								debayerParameters,
 								colorProcParameters_aux,
 								rgbParameters,

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

@@ -144,6 +144,10 @@ public class XyzAtr {
 		atr[2] = d[2];
 	}
 	
+	public double [][] getErsXYZATR_dt() {
+		return new double [][] {ers_xyz_dt, ers_atr_dt};
+	}
+	
 	public double [] getErsXYZ_dt() {
 		return ers_xyz_dt;
 	}

src/main/java/com/elphel/imagej/tileprocessor/lwoc/LwirWorld.java

@@ -340,6 +340,13 @@ public class LwirWorld {
 				quadCLTs[ref_index].saveQuadClt(); // to re-load new set of Bayer images to the GPU (do nothing for CPU) and Geometry
 				quadCLTs[ref_index].preExpandCLTQuad3d( // returns ImagePlus, but it already should be saved/shown
 						clt_parameters,
+						// adding these parameters for more flexibility in accuracy/speed
+						clt_parameters.gr_max_clust_radius, //int      gr_max_clust_radius,			  
+						clt_parameters.disp_scan_start,     // double  disp_scan_start,
+						clt_parameters.disp_scan_step,      // double  disp_scan_step,
+						clt_parameters.disp_scan_count,     // double  disp_scan_count,
+						false,                              // boolean no_bg_generate,			  
+						false,                              // boolean no_lma,			  
 						debayerParameters,
 						colorProcParameters,
 						rgbParameters,

src/main/java/com/elphel/imagej/tileprocessor/sfm/SfmCorr.java

 package com.elphel.imagej.tileprocessor.sfm;
 
 class SfmCorr{
-	double sfm_gain;
-	double [] corr_ind = null;  // {disparity, strength}
+	double sfm_gain =     0;
+	double [] corr_ind =  null;  // {disparity, strength}
 	double [] corr_neib = null; // {disparity, strength}
 }

src/main/java/com/elphel/imagej/tileprocessor/sfm/StructureFromMotion.java

@@ -643,7 +643,7 @@ public class StructureFromMotion {
         				tilesY,
         				true,
         				"sfm_corr-"+scene_pairs[num_pairs-1][1].getImageName()+"-"+
-        						scene_pairs[0][0].getImageName()+"-"+num_pairs,
+        						scene_pairs[0][0].getImageName()+"-"+num_pairs+"-"+String.format("%02d", ntry),
         						dbg_titles);
         	}
         	
@@ -862,6 +862,9 @@ public class StructureFromMotion {
 			threads[ithread] = new Thread() {
 				public void run() {
 					for (int nTile = ai.getAndIncrement(); nTile < accum_PD.length; nTile = ai.getAndIncrement()) {
+//						if (nTile==4600) {
+//							System.out.println("getSfmCorr(): nTile="+nTile);
+//						}
 						if ((disp_corr[0][nTile] != null) || (use_neibs && (disp_corr[1][nTile] != null) )) {
 							sfmCorr[nTile] = new SfmCorr();
 							double [] sxy = dpXYddisp_avg[nTile];
@@ -1007,19 +1010,47 @@ public class StructureFromMotion {
 						debugLevel);                // int            debug_level)
 			}
 		}
+		// remove tiles from both that are null for any
+		final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
+		final AtomicInteger ai = new AtomicInteger(0);
+		for (int ithread = 0; ithread < threads.length; ithread++) {
+			threads[ithread] = new Thread() {
+				public void run() {
+					for	(int nTile = ai.getAndIncrement(); nTile < ref_pXpYD.length; nTile = ai.getAndIncrement())
+						if (ref_pXpYD[nTile] != null){
+							boolean undefined = false;
+							for (int nscene = 0; nscene < scenes.length; nscene++) {
+								if (scenes_pXpYD[nscene][nTile] == null) {
+									undefined=true;
+									break;
+								}
+							}
+							if (undefined) {
+								for (int nscene = 0; nscene < scenes.length; nscene++) {
+									scenes_pXpYD[nscene][nTile] = null;
+								}								
+							}
+						
+					}
+				}
+			};
+		}		      
+		ImageDtt.startAndJoin(threads);
+
+		
 		tp_tasks_ref = // were used only for the FPN removal
 				Interscene.setReferenceGPU (
-				clt_parameters,            // CLTParameters      clt_parameters,			
-				ref_scene,                 // QuadCLT            ref_scene,
-				ref_disparity,             // double []          ref_disparity, // null or alternative reference disparity
-				ref_pXpYD,                 // double [][]        ref_pXpYD,
-				null, // reliable_ref,              // final boolean []   selection, // may be null, if not null do not  process unselected tiles
-				margin,                    // final int          margin,
+				clt_parameters,      // CLTParameters      clt_parameters,			
+				scenes[0],           // ref_scene,                 // QuadCLT            ref_scene,
+				null,                // ref_disparity,             // double []          ref_disparity, // null or alternative reference disparity
+				scenes_pXpYD[0],     // ref_pXpYD,                 // double [][]        ref_pXpYD,
+				null,                // reliable_ref,              // final boolean []   selection, // may be null, if not null do not  process unselected tiles
+				margin,              // final int          margin,
 				// motion blur compensation 
-				mb_tau,                    // double             mb_tau,      // 0.008; // time constant, sec
-				mb_max_gain,               // double             mb_max_gain, // 5.0;   // motion blur maximal gain (if more - move second point more than a pixel
-				mb_vectors_scenes[0], // mb_vectors_ref,            // double [][]        mb_vectors,  // now [2][ntiles];
-				debugLevel);               // int         debugLevel)
+				mb_tau,              // double             mb_tau,      // 0.008; // time constant, sec
+				mb_max_gain,         // double             mb_max_gain, // 5.0;   // motion blur maximal gain (if more - move second point more than a pixel
+				mb_vectors_scenes[0],// mb_vectors_ref,            // double [][]        mb_vectors,  // now [2][ntiles];
+				debugLevel);         // int         debugLevel)
 		// is it needed here?
 		scenes[0].saveQuadClt(); // to re-load new set of Bayer images to the GPU (do nothing for CPU) and Geometry