simultaneously correcting a quad

src/main/java/EyesisCorrectionParameters.java

@@ -1892,6 +1892,7 @@ public class EyesisCorrectionParameters {
   		public double     vignetting_max    = 0.4; // value in vignetting data to correspond to 1x in the kernel
   		public double     vignetting_range  = 5.0; // do not try to correct vignetting less than vignetting_max/vignetting_range
   		public int        kernel_step =       16;  // source kernels step in pixels (have 1 kernel margin on each side)  
+  		public double     disparity  =        0.0; // nominal disparity between side of square cameras (pix) 
   		
   		public CLTParameters(){}
   		public void setProperties(String prefix,Properties properties){
@@ -1919,6 +1920,7 @@ public class EyesisCorrectionParameters {
   			properties.setProperty(prefix+"vignetting_max",   this.vignetting_max+"");
   			properties.setProperty(prefix+"vignetting_range", this.vignetting_range+"");
   			properties.setProperty(prefix+"kernel_step",      this.kernel_step+"");
+  			properties.setProperty(prefix+"disparity",        this.disparity +"");
   			
   		}
   		public void getProperties(String prefix,Properties properties){
@@ -1945,6 +1947,7 @@ public class EyesisCorrectionParameters {
   			if (properties.getProperty(prefix+"vignetting_max")!=null) this.vignetting_max=Double.parseDouble(properties.getProperty(prefix+"vignetting_max"));
   			if (properties.getProperty(prefix+"vignetting_range")!=null) this.vignetting_range=Double.parseDouble(properties.getProperty(prefix+"vignetting_range"));
   			if (properties.getProperty(prefix+"kernel_step")!=null)    this.kernel_step=Integer.parseInt(properties.getProperty(prefix+"kernel_step"));
+  			if (properties.getProperty(prefix+"disparity")!=null)      this.disparity=Double.parseDouble(properties.getProperty(prefix+"disparity"));
   		}
   		
   		public boolean showDialog() {
@@ -1972,6 +1975,7 @@ public class EyesisCorrectionParameters {
   			gd.addNumericField("Value (max) in vignetting data to correspond to 1x in the kernel",        this.vignetting_max,      3);
   			gd.addNumericField("Do not try to correct vignetting smaller than this fraction of max",      this.vignetting_range,  3);
   			gd.addNumericField("Kernel step in pixels (has 1 kernel margin on each side)",                this.kernel_step,            0);
+  			gd.addNumericField("Nominal (rectilinear) disparity between side of square cameras (pix)",    this.disparity,  3);
    			
   			WindowTools.addScrollBars(gd);
   			gd.showDialog();
@@ -2000,6 +2004,7 @@ public class EyesisCorrectionParameters {
   			this.vignetting_max=        gd.getNextNumber();
   			this.vignetting_range=      gd.getNextNumber();
   			this.kernel_step=     (int) gd.getNextNumber();
+  			this.disparity=             gd.getNextNumber();
   			
   			return true;
   		}

src/main/java/EyesisDCT.java

@@ -24,7 +24,6 @@
 
 import java.util.ArrayList;
 import java.util.concurrent.atomic.AtomicInteger;
-
 import ij.CompositeImage;
 import ij.IJ;
 import ij.ImagePlus;
@@ -40,6 +39,7 @@ public class EyesisDCT {
 	public EyesisCorrectionParameters.DCTParameters dctParameters = null;
 	public DCTKernels [] kernels = null;
 	double [][][][][][] clt_kernels = null;
+	GeometryCorrection geometryCorrection = null;
 	public int extra_items = 8; // number of extra items saved with kernels (center offset (partial, full, derivatives)
 	public ImagePlus eyesisKernelImage = null;
 	public long startTime;
@@ -88,7 +88,97 @@ public class EyesisDCT {
 	public boolean CLTKernelsAvailable(){
 		return clt_kernels != null;
 	}
-
+	public boolean geometryCorrectionAvailable(){
+		return geometryCorrection != null;
+	}
+	public boolean initGeometryCorrection(int debugLevel){
+		geometryCorrection = new GeometryCorrection();
+		PixelMapping.SensorData [] sensors =  eyesisCorrections.pixelMapping.sensors;
+		// verify that all sensors have the same distortion parameters
+		int numSensors = sensors.length;
+		for (int i = 1; i < numSensors; i++){
+			if (	(sensors[0].focalLength !=           sensors[i].focalLength) ||
+					(sensors[0].distortionC !=           sensors[i].distortionC) ||
+					(sensors[0].distortionB !=           sensors[i].distortionB) ||
+					(sensors[0].distortionA !=           sensors[i].distortionA) ||
+					(sensors[0].distortionA5 !=          sensors[i].distortionA5) ||
+					(sensors[0].distortionA6 !=          sensors[i].distortionA6) ||
+					(sensors[0].distortionA7 !=          sensors[i].distortionA7) ||
+					(sensors[0].distortionA8 !=          sensors[i].distortionA8) ||
+					(sensors[0].distortionRadius !=      sensors[i].distortionRadius) ||
+					(sensors[0].pixelCorrectionWidth !=  sensors[i].pixelCorrectionWidth) ||
+					(sensors[0].pixelCorrectionHeight != sensors[i].pixelCorrectionHeight) ||
+					(sensors[0].pixelSize !=             sensors[i].pixelSize)){
+				System.out.println("initGeometryCorrection(): All sensors have to have the same distortion model, but channels 0 and "+i+" mismatch");
+				return false;
+			}
+		}
+		// set common distportion parameters
+		geometryCorrection.setDistortion(
+				sensors[0].focalLength,
+				sensors[0].distortionC,
+				sensors[0].distortionB,
+				sensors[0].distortionA,
+				sensors[0].distortionA5,
+				sensors[0].distortionA6,
+				sensors[0].distortionA7,
+				sensors[0].distortionA8,
+				sensors[0].distortionRadius,
+				sensors[0].pixelCorrectionWidth,   // virtual camera center is at (pixelCorrectionWidth/2, pixelCorrectionHeight/2)
+				sensors[0].pixelCorrectionHeight,
+				sensors[0].pixelSize);
+		// set other/individual sensor parameters
+		for (int i = 1; i < numSensors; i++){
+			if (	(sensors[0].theta !=                 sensors[i].theta) || // elevation
+					(sensors[0].heading !=               sensors[i].heading)){
+				System.out.println("initGeometryCorrection(): All sensors have to have the same elevation and heading, but channels 0 and "+i+" mismatch");
+				return false;
+			}
+		}
+		double []   forward = new double[numSensors];
+		double []   right =   new double[numSensors];
+		double []   height =  new double[numSensors];
+		double []   roll =    new double[numSensors];
+		double [][] pXY0 =    new double[numSensors][2];
+		for (int i = 0; i < numSensors; i++){
+			forward[i] = sensors[i].forward;
+			right[i] =   sensors[i].right;
+			height[i] =  sensors[i].height;
+			roll[i] =    sensors[i].psi;
+			pXY0[i][0] = sensors[i].px0;
+			pXY0[i][1] = sensors[i].py0;
+		}
+		geometryCorrection.setSensors(
+				numSensors,
+				sensors[0].theta,
+				sensors[0].heading,
+				forward,
+				right,
+				height,
+				roll, 
+				pXY0);
+		geometryCorrection.planeProjectLenses(); // project all lenses to the common plane 
+		
+		// calcualte reverse distortion as a table to be linear intr4epolated
+		geometryCorrection.calcReverseDistortionTable();
+		
+		if (numSensors == 4){
+			geometryCorrection.adustSquare();
+			System.out.println("Adjusted camera to orient X Y along the sides of a square");
+		} else {
+			System.out.println("============= Cannot adustSquare() as it requires exactly 4 sensors, "+numSensors+" provided ==========");
+			return false;
+		}
+		// Print parameters
+		if (debugLevel > 0){
+			geometryCorrection.listGeometryCorrection(debugLevel > 1);
+		}
+		
+//listGeometryCorrection		
+		return true;
+	}
+	
+//GeometryCorrection
 
 	  public DCTKernels calculateDCTKernel (
 			  final ImageStack kernelStack,  // first stack with 3 colors/slices convolution kernels
@@ -1123,9 +1213,11 @@ public class EyesisDCT {
 	  {
 		  kernels = null;
 	  }
-	  public void resetCLTKernels()
+	  public void resetCLTKernels() // and geometry corection too
 	  {
 		  clt_kernels = null;
+		  geometryCorrection=null;
+
 	  }
 
 	  public boolean readDCTKernels(
@@ -3879,5 +3971,625 @@ public class EyesisDCT {
 	  }
 	  
 
+	  
+	  
+	  
+	  
+	  
+	  public void processCLTQuads(
+			  EyesisCorrectionParameters.CLTParameters           clt_parameters,
+			  EyesisCorrectionParameters.DebayerParameters     debayerParameters,
+			  EyesisCorrectionParameters.NonlinParameters       nonlinParameters,
+			  EyesisCorrectionParameters.ColorProcParameters colorProcParameters,
+			  CorrectionColorProc.ColorGainsParameters     channelGainParameters,
+			  EyesisCorrectionParameters.RGBParameters             rgbParameters,
+			  EyesisCorrectionParameters.EquirectangularParameters equirectangularParameters,
+			  int          convolveFFTSize, // 128 - fft size, kernel size should be size/2
+			  final int          threadsMax,  // maximal number of threads to launch                         
+			  final boolean    updateStatus,
+			  final int        debugLevel)
+	  {
+		  this.startTime=System.nanoTime();
+		  String [] sourceFiles=correctionsParameters.getSourcePaths();
+		  boolean [] enabledFiles=new boolean[sourceFiles.length];
+		  for (int i=0;i<enabledFiles.length;i++) enabledFiles[i]=false;
+		  int numFilesToProcess=0;
+		  int numImagesToProcess=0;
+		  for (int nFile=0;nFile<enabledFiles.length;nFile++){
+			  if ((sourceFiles[nFile]!=null) && (sourceFiles[nFile].length()>1)) {
+				  int [] channels={correctionsParameters.getChannelFromSourceTiff(sourceFiles[nFile])};
+				  if (correctionsParameters.isJP4()){
+					  int subCamera= channels[0]- correctionsParameters.firstSubCamera; // to match those in the sensor files
+					  // removeUnusedSensorData should be off!?
+					  channels=this.eyesisCorrections.pixelMapping.channelsForSubCamera(subCamera);
+				  }
+				  if (channels!=null){
+					  for (int i=0;i<channels.length;i++) if (eyesisCorrections.isChannelEnabled(channels[i])){
+						  if (!enabledFiles[nFile]) numFilesToProcess++;
+						  enabledFiles[nFile]=true;
+						  numImagesToProcess++;
+					  }
+				  }
+			  }
+		  }
+		  if (numFilesToProcess==0){
+			  System.out.println("No files to process (of "+sourceFiles.length+")");
+			  return;
+		  } else {
+			  if (debugLevel>0) System.out.println(numFilesToProcess+ " files to process (of "+sourceFiles.length+"), "+numImagesToProcess+" images to process");
+		  }
+		  double [] referenceExposures=eyesisCorrections.calcReferenceExposures(debugLevel); // multiply each image by this and divide by individual (if not NaN)
+		  int [][] fileIndices=new int [numImagesToProcess][2]; // file index, channel number
+		  int index=0;
+		  for (int nFile=0;nFile<enabledFiles.length;nFile++){
+			  if ((sourceFiles[nFile]!=null) && (sourceFiles[nFile].length()>1)) {
+				  int [] channels={correctionsParameters.getChannelFromSourceTiff(sourceFiles[nFile])};
+				  if (correctionsParameters.isJP4()){
+					  int subCamera= channels[0]- correctionsParameters.firstSubCamera; // to match those in the sensor files
+					  channels=eyesisCorrections.pixelMapping.channelsForSubCamera(subCamera);
+				  }
+				  if (channels!=null){
+					  for (int i=0;i<channels.length;i++) if (eyesisCorrections.isChannelEnabled(channels[i])){
+						  fileIndices[index  ][0]=nFile;
+						  fileIndices[index++][1]=channels[i];
+					  }
+				  }
+			  }
+		  }
+		  ArrayList<String> setNames = new ArrayList<String>();
+		  ArrayList<ArrayList<Integer>> setFiles = new ArrayList<ArrayList<Integer>>();
+
+		  for (int iImage=0;iImage<fileIndices.length;iImage++){
+			  int nFile=fileIndices[iImage][0];
+			  String setName = correctionsParameters.getNameFromSourceTiff(sourceFiles[nFile]);
+			  if (!setNames.contains(setName)) {
+				  setNames.add(setName);
+				  setFiles.add(new ArrayList<Integer>());
+			  }
+			  setFiles.get(setNames.indexOf(setName)).add(new Integer(nFile));
+		  }		  
+		  int iImage = 0;
+		  for (int nSet = 0; nSet < setNames.size(); nSet++){
+			  int maxChn = 0;
+			  for (int i = 0; i < setFiles.get(nSet).size(); i++){
+				  int chn = fileIndices[setFiles.get(nSet).get(i)][1];
+				  if (chn > maxChn) maxChn = chn;
+			  }
+			  int [] channelFiles = new int[maxChn+1];
+			  for (int i =0; i < channelFiles.length; i++) channelFiles[i] = -1;
+			  for (int i = 0; i < setFiles.get(nSet).size(); i++){
+				  channelFiles[fileIndices[setFiles.get(nSet).get(i)][1]] = setFiles.get(nSet).get(i);
+			  }
+			  
+			  ImagePlus [] imp_srcs = new ImagePlus[channelFiles.length];
+			  double [] scaleExposure = new double[channelFiles.length];
+			  for (int srcChannel=0; srcChannel<channelFiles.length; srcChannel++){
+				  int nFile=channelFiles[srcChannel];
+				  imp_srcs[srcChannel]=null;
+				  if (nFile >=0){
+					  if (correctionsParameters.isJP4()){
+						  int subchannel=eyesisCorrections.pixelMapping.getSubChannel(srcChannel);
+						  if (this.correctionsParameters.swapSubchannels01) {
+							  switch (subchannel){
+							  case 0: subchannel=1; break;
+							  case 1: subchannel=0; break;
+							  }
+						  }
+						  if (debugLevel>0) System.out.println("Processing set " + setNames.get(nSet)+" channel "+srcChannel+" - subchannel "+subchannel+" of "+sourceFiles[nFile]);
+						  ImagePlus imp_composite=eyesisCorrections.JP4_INSTANCE.open(
+								  "", // path,
+								  sourceFiles[nFile],
+								  "",  //arg - not used in JP46 reader
+								  true, // un-apply camera color gains
+								  null, // new window
+								  false); // do not show
+						  imp_srcs[srcChannel]=eyesisCorrections.JP4_INSTANCE.demuxImage(imp_composite, subchannel);
+						  if (imp_srcs[srcChannel] == null) imp_srcs[srcChannel] = imp_composite; // not a composite image
+						  // do we need to add any properties?				  
+					  } else { 
+						  imp_srcs[srcChannel]=new ImagePlus(sourceFiles[nFile]);
+						  //					  (new JP46_Reader_camera(false)).decodeProperiesFromInfo(imp_src); // decode existent properties from info
+						  eyesisCorrections.JP4_INSTANCE.decodeProperiesFromInfo(imp_srcs[srcChannel]); // decode existent properties from info
+						  if (debugLevel>0) System.out.println("Processing "+sourceFiles[nFile]);
+					  }
+					  scaleExposure[srcChannel] = 1.0;
+					  if (!Double.isNaN(referenceExposures[nFile]) && (imp_srcs[srcChannel].getProperty("EXPOSURE")!=null)){
+						  scaleExposure[srcChannel] = referenceExposures[nFile]/Double.parseDouble((String) imp_srcs[srcChannel].getProperty("EXPOSURE"));
+						  if (debugLevel>0) System.out.println("Will scale intensity (to compensate for exposure) by  "+scaleExposure);
+					  }
+					  imp_srcs[srcChannel].setProperty("name",    correctionsParameters.getNameFromSourceTiff(sourceFiles[nFile]));
+					  imp_srcs[srcChannel].setProperty("channel", srcChannel); // it may already have channel
+					  imp_srcs[srcChannel].setProperty("path",    sourceFiles[nFile]); // it may already have channel
+
+					  if (this.correctionsParameters.pixelDefects && (eyesisCorrections.defectsXY!=null)&& (eyesisCorrections.defectsXY[srcChannel]!=null)){
+						  // apply pixel correction
+						  int numApplied=	eyesisCorrections.correctDefects(
+								  imp_srcs[srcChannel],
+								  srcChannel,
+								  debugLevel);
+						  if ((debugLevel>0) && (numApplied>0)) { // reduce verbosity after verified defect correction works
+							  System.out.println("Corrected "+numApplied+" pixels in "+sourceFiles[nFile]);
+						  }
+					  }
+
+					  if (this.correctionsParameters.vignetting){
+						  if ((eyesisCorrections.channelVignettingCorrection==null) || (srcChannel<0) || (srcChannel>=eyesisCorrections.channelVignettingCorrection.length) || (eyesisCorrections.channelVignettingCorrection[srcChannel]==null)){
+							  System.out.println("No vignetting data for channel "+srcChannel);
+							  return;
+						  }
+						  float [] pixels=(float []) imp_srcs[srcChannel].getProcessor().getPixels();
+						  if (pixels.length!=eyesisCorrections.channelVignettingCorrection[srcChannel].length){
+							  System.out.println("Vignetting data for channel "+srcChannel+" has "+eyesisCorrections.channelVignettingCorrection[srcChannel].length+" pixels, image "+sourceFiles[nFile]+" has "+pixels.length);
+							  return;
+						  }
+						  // TODO: Move to do it once:
+						  double min_non_zero = 0.0;
+						  for (int i=0;i<pixels.length;i++){
+							  double d = eyesisCorrections.channelVignettingCorrection[srcChannel][i];
+							  if ((d > 0.0) && ((min_non_zero == 0) || (min_non_zero > d))){
+								  min_non_zero = d;
+							  }
+						  }
+						  double max_vign_corr = clt_parameters.vignetting_range*min_non_zero;
+
+						  System.out.println("Vignetting data: channel="+srcChannel+", min = "+min_non_zero);
+						  for (int i=0;i<pixels.length;i++){
+							  double d = eyesisCorrections.channelVignettingCorrection[srcChannel][i];
+							  if (d > max_vign_corr) d = max_vign_corr;
+							  pixels[i]*=d;
+						  }
+						  // Scale here, combine with vignetting later?
+						  int width =  imp_srcs[srcChannel].getWidth();
+						  int height = imp_srcs[srcChannel].getHeight();
+						  for (int y = 0; y < height-1; y+=2){
+							  for (int x = 0; x < width-1; x+=2){
+								  pixels[y*width+x        ] *= clt_parameters.scale_g;
+								  pixels[y*width+x+width+1] *= clt_parameters.scale_g;
+								  pixels[y*width+x      +1] *= clt_parameters.scale_r;
+								  pixels[y*width+x+width  ] *= clt_parameters.scale_b;
+							  }
+						  }
+
+					  } else { // assuming GR/BG pattern
+						  System.out.println("Applying fixed color gain correction parameters: Gr="+
+								  clt_parameters.novignetting_r+", Gg="+clt_parameters.novignetting_g+", Gb="+clt_parameters.novignetting_b);
+						  float [] pixels=(float []) imp_srcs[srcChannel].getProcessor().getPixels();
+						  int width =  imp_srcs[srcChannel].getWidth();
+						  int height = imp_srcs[srcChannel].getHeight();
+						  double kr = clt_parameters.scale_r/clt_parameters.novignetting_r;
+						  double kg = clt_parameters.scale_g/clt_parameters.novignetting_g;
+						  double kb = clt_parameters.scale_b/clt_parameters.novignetting_b;
+						  for (int y = 0; y < height-1; y+=2){
+							  for (int x = 0; x < width-1; x+=2){
+								  pixels[y*width+x        ] *= kg;
+								  pixels[y*width+x+width+1] *= kg;
+								  pixels[y*width+x      +1] *= kr;
+								  pixels[y*width+x+width  ] *= kb;
+							  }
+						  }
+					  }
+				  }
+			  }
+			  // may need to equalize gains between channels
+			  if (clt_parameters.gains_equalize){
+				  double [][] avr_pix = new double [channelFiles.length][3];
+				  double [] avr_RGB = {0.0,0.0,0.0};
+				  int numChn = 0;
+				  for (int srcChannel=0; srcChannel < channelFiles.length; srcChannel++){
+					  int nFile=channelFiles[srcChannel];
+					  if (nFile >=0){
+						  for (int i = 0; i < avr_pix[srcChannel].length; i++) avr_pix[srcChannel][i] = 0;
+						  float [] pixels=(float []) imp_srcs[srcChannel].getProcessor().getPixels();
+						  int width =  imp_srcs[srcChannel].getWidth();
+						  int height = imp_srcs[srcChannel].getHeight();
+						  for (int y = 0; y < height-1; y+=2){
+							  for (int x = 0; x < width-1; x+=2){
+								  avr_pix[srcChannel][0] += pixels[y*width+x      +1]; 
+								  avr_pix[srcChannel][2] += pixels[y*width+x+width  ];
+								  avr_pix[srcChannel][1] += pixels[y*width+x        ];
+								  avr_pix[srcChannel][1] += pixels[y*width+x+width+1];
+							  }
+						  }
+						  avr_pix[srcChannel][0] /= 0.25*width*height;
+						  avr_pix[srcChannel][1] /= 0.5 *width*height;
+						  avr_pix[srcChannel][2] /= 0.25*width*height;
+						  for (int j=0; j < avr_RGB.length; j++) avr_RGB[j] += avr_pix[srcChannel][j];
+						  numChn++;
+						  if (debugLevel>-1) {
+							  System.out.println("processCLTSets(): set "+ setNames.get(nSet) + " channel "+srcChannel+
+									  " R"+avr_pix[srcChannel][0]+" G"+avr_pix[srcChannel][1]+" B"+avr_pix[srcChannel][2]);
+						  }
+						  
+					  }
+				  }
+				  for (int j=0; j < avr_RGB.length; j++) avr_RGB[j] /= numChn;
+				  if (debugLevel>-1) {
+					  System.out.println("processCLTSets(): set "+ setNames.get(nSet) + "average color values: "+
+							  " R="+avr_RGB[0]+" G=" + avr_RGB[1]+" B=" + avr_RGB[2]);
+				  }
+				  for (int srcChannel=0; srcChannel < channelFiles.length; srcChannel++){
+					  int nFile=channelFiles[srcChannel];
+					  if (nFile >=0){
+						  double [] scales = new double [avr_RGB.length];
+						  for (int j=0;j < scales.length; j++){
+							  scales[j] = avr_RGB[j]/avr_pix[srcChannel][j];
+						  }
+						  float [] pixels=(float []) imp_srcs[srcChannel].getProcessor().getPixels();
+						  int width =  imp_srcs[srcChannel].getWidth();
+						  int height = imp_srcs[srcChannel].getHeight();
+						  for (int y = 0; y < height-1; y+=2){
+							  for (int x = 0; x < width-1; x+=2){
+								  pixels[y*width+x        ] *= scales[1];
+								  pixels[y*width+x+width+1] *= scales[1];
+								  pixels[y*width+x      +1] *= scales[0];
+								  pixels[y*width+x+width  ] *= scales[2];
+							  }
+						  }
+					  }
+				  }
+			  }
+
+			  for (int srcChannel=0; srcChannel<channelFiles.length; srcChannel++){
+				  int nFile=channelFiles[srcChannel];
+				  if (nFile >=0){
+					  processCLTQuad( // returns ImagePlus, but it already should be saved/shown
+							  imp_srcs, // [srcChannel], // should have properties "name"(base for saving results), "channel","path"
+							  clt_parameters,
+							  debayerParameters,
+							  nonlinParameters,
+							  colorProcParameters,
+							  channelGainParameters,
+							  rgbParameters,
+							  convolveFFTSize, // 128 - fft size, kernel size should be size/2
+							  scaleExposure[srcChannel],
+							  threadsMax,  // maximal number of threads to launch                         
+							  updateStatus,
+							  debugLevel);
+					  // warp result (add support for different color modes)
+					  if (this.correctionsParameters.equirectangular){
+						  if (equirectangularParameters.clearFullMap) eyesisCorrections.pixelMapping.deleteEquirectangularMapFull(srcChannel); // save memory? //removeUnusedSensorData - no, use equirectangular specific settings
+						  if (equirectangularParameters.clearAllMaps) eyesisCorrections.pixelMapping.deleteEquirectangularMapAll(srcChannel); // save memory? //removeUnusedSensorData - no, use equirectangular specific settings
+					  }
+					  //pixelMapping
+					  Runtime.getRuntime().gc();
+					  if (debugLevel >-1) System.out.println("Processing image "+(iImage+1)+" (of "+fileIndices.length+") finished at "+
+							  IJ.d2s(0.000000001*(System.nanoTime()-this.startTime),3)+" sec, --- Free memory="+Runtime.getRuntime().freeMemory()+" (of "+Runtime.getRuntime().totalMemory()+")");
+					  if (eyesisCorrections.stopRequested.get()>0) {
+						  System.out.println("User requested stop");
+						  return;
+					  }
+					  iImage++;
+				  }
+			  }
+		  }
+		  System.out.println("Processing "+fileIndices.length+" files finished at "+
+				  IJ.d2s(0.000000001*(System.nanoTime()-this.startTime),3)+" sec, --- Free memory="+Runtime.getRuntime().freeMemory()+" (of "+Runtime.getRuntime().totalMemory()+")");
+	  }		
+	  
+	  public ImagePlus [] processCLTQuad(
+			  ImagePlus [] imp_quad, // should have properties "name"(base for saving results), "channel","path"
+			  EyesisCorrectionParameters.CLTParameters           clt_parameters,
+			  EyesisCorrectionParameters.DebayerParameters     debayerParameters,
+			  EyesisCorrectionParameters.NonlinParameters       nonlinParameters,
+			  EyesisCorrectionParameters.ColorProcParameters colorProcParameters,
+			  CorrectionColorProc.ColorGainsParameters     channelGainParameters,
+			  EyesisCorrectionParameters.RGBParameters             rgbParameters,
+			  int          convolveFFTSize, // 128 - fft size, kernel size should be size/2
+			  double 		     scaleExposure,
+			  final int        threadsMax,  // maximal number of threads to launch                         
+			  final boolean    updateStatus,
+			  final int        debugLevel){
+		  boolean advanced=this.correctionsParameters.zcorrect || this.correctionsParameters.equirectangular;
+//		  boolean crop=      advanced? true: this.correctionsParameters.crop; 
+		  boolean rotate=    advanced? false: this.correctionsParameters.rotate; 
+		  double JPEG_scale= advanced? 1.0: this.correctionsParameters.JPEG_scale;
+		  boolean toRGB=     advanced? true: this.correctionsParameters.toRGB; 
+		  showDoubleFloatArrays sdfa_instance = new showDoubleFloatArrays(); // just for debugging?
+
+		  // may use this.StartTime to report intermediate steps execution times
+		  String name=(String) imp_quad[0].getProperty("name");
+		  //		int channel= Integer.parseInt((String) imp_src.getProperty("channel"));
+		  int channel= (Integer) imp_quad[0].getProperty("channel");
+		  String path= (String) imp_quad[0].getProperty("path");
+		  		  
+//		  String title=name+"-"+String.format("%02d", channel);
+//		  ImagePlus result=imp_src;
+		  ImagePlus [] results = new ImagePlus[imp_quad.length];
+		  for (int i = 0; i < results.length; i++) {
+			  results[i] = imp_quad[i];
+			  results[i].setTitle(results[i].getTitle()+"RAW");			  
+		  }
+		  if (debugLevel>1) System.out.println("processing: "+path);
+//		  result.setTitle(title+"RAW");
+		  double [][][] double_stacks = new double [imp_quad.length][][];
+		  for (int i = 0; i < double_stacks.length; i++){
+			  double_stacks[i] = eyesisCorrections.bayerToDoubleStack(
+					  imp_quad[i], // source Bayer image, linearized, 32-bit (float))
+					  null); // no margins, no oversample
+			  
+		  }
+		  
+		  String [] rbg_titles = {"Red", "Blue", "Green"};
+		  ImageStack stack;
+		  // =================
+		  ImageDtt image_dtt = new ImageDtt();
+		  for (int i = 0; i < double_stacks.length; i++){
+			  for (int j =0 ; j < double_stacks[i][0].length; j++){
+				  double_stacks[i][2][j]*=0.5; // Scale green 0.5 to compensate more pixels than R,B
+			  }
+		  }
+		  double [][][][][][] clt_data = image_dtt.clt_aberrations_quad( 
+				  clt_parameters.disparity,     // final double            disparity,
+				  double_stacks,                // final double [][][]      imade_data, // first index - number of image in a quad
+				  imp_quad[0].getWidth(),       //	final int               width,
+				  geometryCorrection,           // final GeometryCorrection  geometryCorrection,
+				  clt_kernels,                  // final double [][][][][][] clt_kernels, // [channel_in_quad][color][tileY][tileX][band][pixel] , size should match image (have 1 tile around)
+				  clt_parameters.kernel_step,
+				  clt_parameters.transform_size,
+				  clt_parameters.clt_window,
+				  clt_parameters.shift_x,       // final int               shiftX, // shift image horizontally (positive - right) - just for testing
+				  clt_parameters.shift_y,       // final int               shiftY, // shift image vertically (positive - down)
+				  clt_parameters.tileX,         // final int               debug_tileX,
+				  clt_parameters.tileY,         // final int               debug_tileY,
+				  (clt_parameters.dbg_mode & 64) != 0, // no fract shift
+				  (clt_parameters.dbg_mode & 128) != 0, // no convolve
+				  (clt_parameters.dbg_mode & 256) != 0, // transpose convolve
+				  threadsMax,
+				  debugLevel);
+		  //					  updateStatus);
+
+
+		  System.out.println("clt_data.length="+clt_data.length+" clt_data[0].length="+clt_data[0].length
+				  +" clt_data[0][0].length="+clt_data[0][0].length+" clt_data[0][0][0].length="+
+				  clt_data[0][0][0].length+" clt_data[0][0][0][0].length="+clt_data[0][0][0][0].length+
+				  " clt_data[0][0][0][0][0].length="+clt_data[0][0][0][0][0].length);
+		  
+		  for (int iQuad = 0; iQuad <clt_data.length; iQuad++){
+			  
+			  String title=name+"-"+String.format("%02d", iQuad);
+			  String titleFull=title+"-SPLIT";
+			  
+			  if (clt_parameters.corr_sigma > 0){ // no filter at all
+				  for (int chn = 0; chn < clt_data[iQuad].length; chn++) {
+					  image_dtt.clt_lpf(
+							  clt_parameters.corr_sigma,
+							  clt_data[iQuad][chn],
+							  threadsMax,
+							  debugLevel);
+				  }
+			  }
+
+			  int tilesY = imp_quad[iQuad].getHeight()/clt_parameters.transform_size;
+			  int tilesX = imp_quad[iQuad].getWidth()/clt_parameters.transform_size;
+			  if (debugLevel > 0){
+				  System.out.println("--tilesX="+tilesX);
+				  System.out.println("--tilesY="+tilesY);
+			  }
+			  if (debugLevel > 1){
+				  double [][] clt = new double [clt_data[iQuad].length*4][];
+				  for (int chn = 0; chn < clt_data[iQuad].length; chn++) {
+					  double [][] clt_set = image_dtt.clt_dbg(
+							  clt_data [iQuad][chn],
+							  threadsMax,
+							  debugLevel);
+					  for (int ii = 0; ii < clt_set.length; ii++) clt[chn*4+ii] = clt_set[ii];
+				  }
+
+				  if (debugLevel > 0){
+					  sdfa_instance.showArrays(clt,
+							  tilesX*clt_parameters.transform_size,
+							  tilesY*clt_parameters.transform_size,
+							  true,
+							  results[iQuad].getTitle()+"-CLT");  
+				  }
+			  }
+			  double [][] iclt_data = new double [clt_data[iQuad].length][];
+			  for (int chn=0; chn<iclt_data.length;chn++){
+				  iclt_data[chn] = image_dtt.iclt_2d(
+						  clt_data[iQuad][chn],           // scanline representation of dcd data, organized as dct_size x dct_size tiles  
+						  clt_parameters.transform_size,  // final int
+						  clt_parameters.clt_window,      // window_type
+						  15,                             // clt_parameters.iclt_mask,       //which of 4 to transform back
+						  0,                              // clt_parameters.dbg_mode,        //which of 4 to transform back
+						  threadsMax,
+						  debugLevel);
+
+			  }
+			  if (debugLevel > 0) sdfa_instance.showArrays(
+					  iclt_data,
+					  (tilesX + 0) * clt_parameters.transform_size,
+					  (tilesY + 0) * clt_parameters.transform_size,
+					  true,
+					  results[iQuad].getTitle()+"-rbg_sigma");
+			  /*	  
+				  }
+			  }
+			   */
+			  if (debugLevel > 0) sdfa_instance.showArrays(iclt_data,
+					  (tilesX + 0) * clt_parameters.transform_size,
+					  (tilesY + 0) * clt_parameters.transform_size,
+					  true,
+					  results[iQuad].getTitle()+"-ICLT-RGB");
+
+			  // convert to ImageStack of 3 slices
+			  String [] sliceNames = {"red", "blue", "green"};
+			  stack = sdfa_instance.makeStack(
+					  iclt_data,
+					  (tilesX + 0) * clt_parameters.transform_size,
+					  (tilesY + 0) * clt_parameters.transform_size,
+					  sliceNames); // or use null to get chn-nn slice names
+
+
+
+			  if (debugLevel > -1){
+				  double [] chn_avg = {0.0,0.0,0.0};
+				  float [] pixels;
+				  int width =  stack.getWidth();
+				  int height = stack.getHeight();
+
+				  for (int c = 0; c <3; c++){
+					  pixels = (float[]) stack.getPixels(c+1);
+					  for (int i = 0; i<pixels.length; i++){
+						  chn_avg[c] += pixels[i];
+					  }
+				  }
+				  chn_avg[0] /= width*height;
+				  chn_avg[1] /= width*height;
+				  chn_avg[2] /= width*height;
+				  System.out.println("Processed channels averages: R="+chn_avg[0]+", G="+chn_avg[2]+", B="+chn_avg[1]); 
+			  }
+
+			  if (!this.correctionsParameters.colorProc){
+				  results[iQuad]= new ImagePlus(titleFull, stack);    			  
+				  eyesisCorrections.saveAndShow(
+						  results[iQuad],
+						  this.correctionsParameters);
+				  continue; // return results;
+			  }
+			  if (debugLevel > 1) System.out.println("before colors.1");
+			  //Processing colors - changing stack sequence to r-g-b (was r-b-g)
+			  if (!eyesisCorrections.fixSliceSequence(
+					  stack,
+					  debugLevel)){
+				  if (debugLevel > -1) System.out.println("fixSliceSequence() returned false");
+				  return null;
+			  }
+			  if (debugLevel > 1) System.out.println("before colors.2");
+			  if (debugLevel > 1){
+				  ImagePlus imp_dbg=new ImagePlus(imp_quad[iQuad].getTitle()+"-"+channel+"-preColors",stack);
+				  eyesisCorrections.saveAndShow(
+						  imp_dbg,
+						  this.correctionsParameters);
+			  }
+			  if (debugLevel > 1) System.out.println("before colors.3, scaleExposure="+scaleExposure+" scale = "+(255.0/eyesisCorrections.psfSubpixelShouldBe4/eyesisCorrections.psfSubpixelShouldBe4/scaleExposure));
+			  CorrectionColorProc correctionColorProc=new CorrectionColorProc(eyesisCorrections.stopRequested);
+			  double [][] yPrPb=new double [3][];
+			  //			if (dct_parameters.color_DCT){
+			  // need to get YPbPr - not RGB here				
+			  //			} else {
+			  correctionColorProc.processColorsWeights(stack, // just gamma convert? TODO: Cleanup? Convert directly form the linear YPrPb
+					  //					  255.0/this.psfSubpixelShouldBe4/this.psfSubpixelShouldBe4, //  double scale,     // initial maximal pixel value (16))
+					  //					  255.0/eyesisCorrections.psfSubpixelShouldBe4/eyesisCorrections.psfSubpixelShouldBe4/scaleExposure, //  double scale,     // initial maximal pixel value (16))
+					  //					  255.0/2/2/scaleExposure, //  double scale,     // initial maximal pixel value (16))
+					  255.0/scaleExposure, //  double scale,     // initial maximal pixel value (16))
+					  colorProcParameters,
+					  channelGainParameters,
+					  channel,
+					  null, //correctionDenoise.getDenoiseMask(),
+					  this.correctionsParameters.blueProc,
+					  debugLevel);
+			  if (debugLevel > 1) System.out.println("Processed colors to YPbPr, total number of slices="+stack.getSize());
+			  if (debugLevel > 1) {
+				  ImagePlus imp_dbg=new ImagePlus("procColors",stack);
+				  eyesisCorrections.saveAndShow(
+						  imp_dbg,
+						  this.correctionsParameters);
+			  }
+			  float [] fpixels;
+			  int [] slices_YPrPb = {8,6,7};
+			  yPrPb=new double [3][];
+			  for (int n = 0; n < slices_YPrPb.length; n++){
+				  fpixels = (float[]) stack.getPixels(slices_YPrPb[n]);
+				  yPrPb[n] = new double [fpixels.length];
+				  for (int i = 0; i < fpixels.length; i++) yPrPb[n][i] = fpixels[i];
+			  }
+
+			  if (toRGB) {
+				  System.out.println("correctionColorProc.YPrPbToRGB");
+				  stack =  YPrPbToRGB(yPrPb,
+						  colorProcParameters.kr,        // 0.299;
+						  colorProcParameters.kb,        // 0.114;
+						  stack.getWidth());
+
+				  title=titleFull; // including "-DECONV" or "-COMBO"
+				  titleFull=title+"-RGB-float";
+				  //Trim stack to just first 3 slices
+				  if (debugLevel > 1){ // 2){
+					  ImagePlus imp_dbg=new ImagePlus("YPrPbToRGB",stack);
+					  eyesisCorrections.saveAndShow(
+							  imp_dbg,
+							  this.correctionsParameters);
+				  }
+				  while (stack.getSize() > 3) stack.deleteLastSlice();
+				  if (debugLevel > 1) System.out.println("Trimming color stack");
+			  } else {
+				  title=titleFull; // including "-DECONV" or "-COMBO"
+				  titleFull=title+"-YPrPb"; // including "-DECONV" or "-COMBO"
+				  if (debugLevel > 1) System.out.println("Using full stack, including YPbPr");
+			  }
+
+			  results[iQuad]= new ImagePlus(titleFull, stack);    			  
+			  // rotate the result			  
+			  if (rotate){ // never rotate for equirectangular
+				  stack=eyesisCorrections.rotateStack32CW(stack);
+			  }
+			  if (!toRGB && !this.correctionsParameters.jpeg){ // toRGB set for equirectangular
+				  if (debugLevel > 1) System.out.println("!toRGB && !this.correctionsParameters.jpeg");
+				  eyesisCorrections.saveAndShow(results[iQuad], this.correctionsParameters);
+				  continue; // return result;
+			  } else { // that's not the end result, save if required
+				  if (debugLevel > 1) System.out.println("!toRGB && !this.correctionsParameters.jpeg - else");
+				  eyesisCorrections.saveAndShow(results[iQuad],
+						  eyesisCorrections.correctionsParameters,
+						  eyesisCorrections.correctionsParameters.save32,
+						  false,
+						  eyesisCorrections.correctionsParameters.JPEG_quality); // save, no show
+			  }
+			  // convert to RGB48 (16 bits per color component)
+			  ImagePlus imp_RGB;
+			  stack=eyesisCorrections.convertRGB32toRGB16Stack(
+					  stack,
+					  rgbParameters); 
+
+			  titleFull=title+"-RGB48";
+			  results[iQuad]= new ImagePlus(titleFull, stack);    			  
+			  
+			  //			  ImagePlus imp_RGB24;
+			  results[iQuad].updateAndDraw();
+			  if (debugLevel > 1) System.out.println("result.updateAndDraw(), "+titleFull+"-RGB48");
+
+			  CompositeImage compositeImage=eyesisCorrections.convertToComposite(results[iQuad]);
+
+			  if (!this.correctionsParameters.jpeg && !advanced){ // RGB48 was the end result
+				  if (debugLevel > 1) System.out.println("if (!this.correctionsParameters.jpeg && !advanced)");
+				  eyesisCorrections.saveAndShow(compositeImage, this.correctionsParameters);
+				  continue; // return result;
+			  } else { // that's not the end result, save if required
+				  if (debugLevel > 1) System.out.println("if (!this.correctionsParameters.jpeg && !advanced) - else");
+				  eyesisCorrections.saveAndShow(compositeImage, this.correctionsParameters, this.correctionsParameters.save16, false); // save, no show
+				  //				  eyesisCorrections.saveAndShow(compositeImage, this.correctionsParameters, this.correctionsParameters.save16, true); // save, no show
+			  }
+
+			  imp_RGB=eyesisCorrections.convertRGB48toRGB24(
+					  stack,
+					  title+"-RGB24",
+					  0, 65536, // r range 0->0, 65536->256
+					  0, 65536, // g range
+					  0, 65536);// b range
+			  if (JPEG_scale!=1.0){
+				  ImageProcessor ip=imp_RGB.getProcessor();
+				  ip.setInterpolationMethod(ImageProcessor.BICUBIC);
+				  ip=ip.resize((int)(ip.getWidth()*JPEG_scale),(int) (ip.getHeight()*JPEG_scale));
+				  imp_RGB= new ImagePlus(imp_RGB.getTitle(),ip);
+				  imp_RGB.updateAndDraw();
+			  }
+			  eyesisCorrections.saveAndShow(imp_RGB, this.correctionsParameters);
+		  }
+		  return results;
+	  }
+	  
+
+	  
+	  
+	  
+	  
+	  
+	  
+	  
+	  
+	  
+	  
+	  
+	  
 
 }

src/main/java/Eyesis_Correction.java

@@ -495,6 +495,7 @@ private Panel panel1,
 			addButton("Reset CLT kernels",         panelClt1, color_stop);
 			addButton("CLT process files",         panelClt1, color_process);
 			addButton("CLT process sets",          panelClt1, color_process);
+			addButton("CLT process quads",         panelClt1, color_process);
 			add(panelClt1);
 		}
 		pack();
@@ -4282,6 +4283,15 @@ private Panel panel1,
             			DEBUG_LEVEL);
         	}
         }
+
+        if (!EYESIS_DCT.geometryCorrectionAvailable()){
+        	if (DEBUG_LEVEL > 0){
+        		System.out.println("Calculating geometryCorrection");
+        	}
+        	if (!EYESIS_DCT.initGeometryCorrection(DEBUG_LEVEL+2)){
+        		return;
+        	}
+        }
         
 ///========================================        
         
@@ -4307,6 +4317,100 @@ private Panel panel1,
         }
         return;
         
+    } else if (label.equals("CLT process quads")) {
+    	DEBUG_LEVEL=MASTER_DEBUG_LEVEL;
+    	EYESIS_CORRECTIONS.setDebug(DEBUG_LEVEL);
+        if (EYESIS_DCT == null){
+        	EYESIS_DCT = new  EyesisDCT (
+        			EYESIS_CORRECTIONS,
+        			CORRECTION_PARAMETERS,
+        			DCT_PARAMETERS);
+        	if (DEBUG_LEVEL > 0){
+        		System.out.println("Created new EyesisDCT instance, will need to read CLT kernels");
+        	}
+        }
+    	String configPath=null;
+    	if (EYESIS_CORRECTIONS.correctionsParameters.saveSettings) {
+    		configPath=EYESIS_CORRECTIONS.correctionsParameters.selectResultsDirectory(
+    				true,
+    				true);
+    		if (configPath==null){
+    			String msg="No results directory selected, command aborted";
+    			System.out.println("Warning: "+msg);
+    			IJ.showMessage("Warning",msg);
+    			return;
+    		}
+    		configPath+=Prefs.getFileSeparator()+"autoconfig";
+    		try {
+    			saveTimestampedProperties(
+    					configPath,      // full path or null
+    					null, // use as default directory if path==null 
+    					true,
+    					PROPERTIES);
+
+    		} catch (Exception e){
+    			String msg="Failed to save configuration to "+configPath+", command aborted";
+    			System.out.println("Error: "+msg);
+    			IJ.showMessage("Error",msg);
+    			return;
+    		}
+    	}      
+        
+        EYESIS_CORRECTIONS.initSensorFiles(DEBUG_LEVEL);
+        int numChannels=EYESIS_CORRECTIONS.getNumChannels();
+        NONLIN_PARAMETERS.modifyNumChannels(numChannels);
+        CHANNEL_GAINS_PARAMETERS.modifyNumChannels(numChannels);
+
+        if (!EYESIS_DCT.CLTKernelsAvailable()){
+        	if (DEBUG_LEVEL > 0){
+        		System.out.println("Reading CLT kernels");
+        	}
+            EYESIS_DCT.readCLTKernels(
+            		CLT_PARAMETERS,
+                    THREADS_MAX,
+                    UPDATE_STATUS, // update status info
+            		DEBUG_LEVEL);
+
+            if (DEBUG_LEVEL > 1){
+            	EYESIS_DCT.showCLTKernels(
+            			THREADS_MAX,
+            			UPDATE_STATUS, // update status info
+            			DEBUG_LEVEL);
+        	}
+        }
+
+        if (!EYESIS_DCT.geometryCorrectionAvailable()){
+        	if (DEBUG_LEVEL > 0){
+        		System.out.println("Calculating geometryCorrection");
+        	}
+        	if (!EYESIS_DCT.initGeometryCorrection(DEBUG_LEVEL+2)){
+        		return;
+        	}
+        }
+        
+///========================================        
+        
+        EYESIS_DCT.processCLTQuads(
+        		CLT_PARAMETERS,  // EyesisCorrectionParameters.DCTParameters           dct_parameters,
+        		DEBAYER_PARAMETERS, //EyesisCorrectionParameters.DebayerParameters     debayerParameters,
+        		NONLIN_PARAMETERS, //EyesisCorrectionParameters.NonlinParameters       nonlinParameters,
+        		COLOR_PROC_PARAMETERS, //EyesisCorrectionParameters.ColorProcParameters colorProcParameters,
+        		CHANNEL_GAINS_PARAMETERS, //CorrectionColorProc.ColorGainsParameters     channelGainParameters,
+        		RGB_PARAMETERS, //EyesisCorrectionParameters.RGBParameters             rgbParameters,
+        		EQUIRECTANGULAR_PARAMETERS, // EyesisCorrectionParameters.EquirectangularParameters equirectangularParameters,
+        		CONVOLVE_FFT_SIZE, //int          convolveFFTSize, // 128 - fft size, kernel size should be size/2
+        		THREADS_MAX, //final int          threadsMax,  // maximal number of threads to launch                         
+        		UPDATE_STATUS, //final boolean    updateStatus,
+        		DEBUG_LEVEL); //final int        debugLevel);
+        
+        if (configPath!=null) {
+        	saveTimestampedProperties( // save config again
+        			configPath,      // full path or null
+        			null, // use as default directory if path==null 
+        			true,
+        			PROPERTIES);
+        }
+        return;
         
         
     	

src/main/java/GeometryCorrection.java

@@ -114,6 +114,7 @@ public class GeometryCorrection {
 			double      heading,    // NaN - keep
 			double []   forward,    // null - keep all, NaN - keep individual
 			double []   right,      // null - keep all, NaN - keep individual
+			double []   height,      // null - keep all, NaN - keep individual
 			double []   roll,       // null - keep all, NaN - keep individual
 			double [][] pXY0){      // null - keep all, [] null - keep individual
 
@@ -136,6 +137,14 @@ public class GeometryCorrection {
 			for (int i = 0; i < numSensors; i++) if (!Double.isNaN(right[i]))      this.right[i] = right[i];
 		}
 
+		if (height != null){
+			if (height.length != numSensors){
+				throw new IllegalArgumentException ("height.length ("+height.length+") != numSensors ("+numSensors+")");
+			}
+			if ((this.height == null) || (this.height.length != numSensors)) this.height = new double [numSensors];
+			for (int i = 0; i < numSensors; i++) if (!Double.isNaN(height[i]))      this.height[i] = height[i];
+		}
+
 		if (roll != null){
 			if (roll.length != numSensors){
 				throw new IllegalArgumentException ("roll.length ("+roll.length+") != numSensors ("+numSensors+")");
@@ -244,10 +253,61 @@ public class GeometryCorrection {
 			for (int j = 0; j<2;j++) this.rXY[i][j] = this.XYZ_her[i][j]/this.disparityRadius; 
 		}
 	}
+	
+	public void listGeometryCorrection(boolean showAll){
+		System.out.println("'=== Constant parameters ===");
+		System.out.println("pixelCorrectionWidth =\t"+  pixelCorrectionWidth+"\tpix");
+		System.out.println("pixelCorrectionHeight =\t"+ pixelCorrectionHeight+"\tpix");
+		System.out.println("pixelSize =\t"+             pixelSize+"\tum");
+		System.out.println("distortionRadius =\t"+      distortionRadius+"\tmm");
+		System.out.println("'=== Common input parameters ===");
+		System.out.println("focalLength =\t"+  focalLength + " mm");
+		System.out.println("distortionA8 =\t"+ distortionA8);
+		System.out.println("distortionA7 =\t"+ distortionA7);
+		System.out.println("distortionA6 =\t"+ distortionA6);
+		System.out.println("distortionA5 =\t"+ distortionA5);
+		System.out.println("distortionA =\t"+  distortionA);
+		System.out.println("distortionB =\t"+  distortionB);
+		System.out.println("distortionC =\t"+  distortionC);
+		System.out.println("elevation =\t"+    elevation+"\tdegrees");
+		System.out.println("heading =\t"+      heading+"\tdegrees");
+		System.out.println("numSensors =\t"+   numSensors);
+		System.out.println("'=== Individual input parameters ===");
+		System.out.print  ("forward = ");for (int i = 0; i < numSensors;i++) System.out.print("\t"+forward[i]); System.out.println("\tmm");
+		System.out.print  ("right = ");  for (int i = 0; i < numSensors;i++) System.out.print("\t"+right[i]);   System.out.println("\tmm");
+		System.out.print  ("height = "); for (int i = 0; i < numSensors;i++) System.out.print("\t"+height[i]);  System.out.println("\tmm");
+		System.out.print  ("roll = ");   for (int i = 0; i < numSensors;i++) System.out.print("\t"+roll[i]);      System.out.println("\tdegrees");
+		System.out.print  ("px0 = ");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+pXY0[i][0]);  System.out.println("\tpix");
+		System.out.print  ("py0 = ");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+pXY0[i][1]);  System.out.println("\tpix");
+		
+		System.out.println("'=== Common calculated parameters ===");
+		System.out.println("common_right =\t"+common_right + "\tmm");
+		System.out.println("common_forward =\t"+common_forward + "\tmm");
+		System.out.println("common_height =\t"+common_height + "\tmm");
+		System.out.println("common_roll =\t"+common_roll + "\tdegrees");
+		System.out.println("cameraRadius =\t"+cameraRadius + "\tmm");
+		System.out.println("disparityRadius =\t"+disparityRadius + "\tmm");
+
+		if (showAll){
+			System.out.println("'=== Intermediate data: coordinates corrected for common elevation and heading ===");
+			System.out.print  ("X_he =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+XYZ_he[i][0]);  System.out.println("\tmm");
+			System.out.print  ("Y_he =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+XYZ_he[i][1]);  System.out.println("\tmm");
+			System.out.print  ("Z_he =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+XYZ_he[i][2]);  System.out.println("\tmm");
+			System.out.println("'=== Intermediate data: coordinates corrected for common elevation, heading and roll ===");
+			System.out.print  ("X_her =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+XYZ_her[i][0]);  System.out.println("\tmm");
+			System.out.print  ("Y_her =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+XYZ_her[i][1]);  System.out.println("\tmm");
+			System.out.print  ("Z_her =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+XYZ_her[i][2]);  System.out.println("\tmm");
+		}
+		
+		System.out.println("'=== Individual calculated parameters ===");
+		System.out.print  ("residual_roll = ");   for (int i = 0; i < numSensors;i++) System.out.print("\t"+(roll[i]-common_roll));System.out.println("\tdegrees");
+		System.out.print  ("X_rel =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+rXY[i][0]);  System.out.println("\trelative to disparityRadius");
+		System.out.print  ("Y_rel =");    for (int i = 0; i < numSensors;i++) System.out.print("\t"+rXY[i][1]);  System.out.println("\trelative to disparityRadius");
+	}
 
 	// return distance from disparity (in pixel units) for the current camera geometry
 	public double getZFromDisparity(double disparity){
-		return disparity * this.focalLength * this.pixelSize / this.disparityRadius;
+		return disparity * this.focalLength * 0.001*this.pixelSize / this.disparityRadius;
 	}
 
 	/*
@@ -263,8 +323,8 @@ public class GeometryCorrection {
 		double [][] pXY = new double [numSensors][2];
 		double pXcd = px - 0.5 * this.pixelCorrectionWidth;
 		double pYcd = py - 0.5 * this.pixelCorrectionHeight;
-		double rD = Math.sqrt(pXcd*pXcd + pYcd*pYcd); // distorted radius in a virtual center camera
-		double rND2R=getRByRDist(rD/this.distortionRadius, (debugLevel > 1));
+		double rD = Math.sqrt(pXcd*pXcd + pYcd*pYcd)*0.001*this.pixelSize; // distorted radius in a virtual center camera
+		double rND2R=getRByRDist(rD/this.distortionRadius, (debugLevel > -1));
 		double pXc = pXcd * rND2R; // non-distorted coordinates relative to the (0.5 * this.pixelCorrectionWidth, 0.5 * this.pixelCorrectionHeight)
 		double pYc = pYcd * rND2R; // in pixels
 		double [] a={this.distortionC,this.distortionB,this.distortionA,this.distortionA5,this.distortionA6,this.distortionA7,this.distortionA8};
@@ -275,7 +335,7 @@ public class GeometryCorrection {
 			// calculate back to distorted
 			double rNDi = Math.sqrt(pXci*pXci + pYci*pYci); // in pixels
 			//		Rdist/R=A8*R^7+A7*R^6+A6*R^5+A5*R^4+A*R^3+B*R^2+C*R+(1-A6-A7-A6-A5-A-B-C)");
-			double ri = rNDi* this.pixelSize / this.distortionRadius; // relative to distortion radius
+			double ri = rNDi* 0.001 * this.pixelSize / this.distortionRadius; // relative to distortion radius
 			//    		double rD2rND = (1.0 - distortionA8 - distortionA7 - distortionA6 - distortionA5 - distortionA - distortionB - distortionC);
 			double rD2rND = 1.0;
 			double rri = 1.0;
@@ -286,8 +346,8 @@ public class GeometryCorrection {
 			double pXid = pXci * rD2rND;  
 			double pYid = pYci * rD2rND;
 			// individual rotate (check sign)
-			double c_roll = Math.cos(this.roll[i] * Math.PI/180.0);
-			double s_roll = Math.sin(this.roll[i] * Math.PI/180.0);
+			double c_roll = Math.cos((this.roll[i] - this.common_roll) * Math.PI/180.0);
+			double s_roll = Math.sin((this.roll[i] - this.common_roll) * Math.PI/180.0);
 			pXY[i][0] =  c_roll *  pXid + s_roll* pYid + this.pXY0[i][0];
 			pXY[i][1] = -s_roll *  pXid + c_roll* pYid + this.pXY0[i][1];
 		}

src/main/java/ImageDtt.java

@@ -757,6 +757,141 @@ public class ImageDtt {
 	}
 
 
+	public double [][][][][][] clt_aberrations_quad(
+			final double              disparity,
+			final double [][][]       imade_data, // first index - number of image in a quad
+			final int                 width,
+			final GeometryCorrection  geometryCorrection,
+			final double [][][][][][] clt_kernels, // [channel_in_quad][color][tileY][tileX][band][pixel] , size should match image (have 1 tile around)
+			final int                 kernel_step,
+			final int                 transform_size,
+			final int                 window_type,
+			final double              shiftX, // shift image horizontally (positive - right) - just for testing
+			final double              shiftY, // shift image vertically (positive - down)
+			final int                 debug_tileX,
+			final int                 debug_tileY,
+			final boolean             no_fract_shift,
+			final boolean             no_deconvolution,
+			final boolean             transpose,
+			final int                 threadsMax,  // maximal number of threads to launch                         
+			final int                 globalDebugLevel)
+	{
+		final int quad = 4;
+		final int nChn = imade_data[0].length;
+		final int height=imade_data[0][0].length/width;
+		final int tilesX=width/transform_size;
+		final int tilesY=height/transform_size;
+		final int nTilesInChn=tilesX*tilesY; 
+		final int nTiles=tilesX*tilesY*nChn; 
+		final double [][][][][][] clt_data = new double[quad][nChn][tilesY][tilesX][4][];
+		final Thread[] threads = newThreadArray(threadsMax);
+		final AtomicInteger ai = new AtomicInteger(0);
+		
+		double [][] centersXY_dbg = geometryCorrection.getPortsCoordinates(
+				0.0, // centerX,
+				0.0, // centerY,
+				disparity);
+		
+		
+		
+		
+		if (globalDebugLevel > 0) {
+			System.out.println("clt_aberrations(): width="+width+" height="+height+" transform_size="+transform_size+
+					" debug_tileX="+debug_tileX+" debug_tileY="+debug_tileY+" globalDebugLevel="+globalDebugLevel);
+		}
+		for (int ithread = 0; ithread < threads.length; ithread++) {
+			threads[ithread] = new Thread() {
+				public void run() {
+					DttRad2 dtt = new DttRad2(transform_size);
+					dtt.set_window(window_type);
+					int tileY,tileX, chn;
+					//						showDoubleFloatArrays sdfa_instance = new showDoubleFloatArrays(); // just for debugging?
+					double centerX; // center of aberration-corrected (common model) tile, X
+					double centerY; //
+					double [][] fract_shiftsXY = new double[quad][];
+
+					for (int nTile = ai.getAndIncrement(); nTile < nTiles; nTile = ai.getAndIncrement()) {
+						// TODO: make all color channels to be processed here (atomically)
+						chn=nTile/nTilesInChn;
+						tileY =(nTile % nTilesInChn)/tilesX;
+						tileX = nTile % tilesX;
+						centerX = tileX * transform_size + transform_size/2 - shiftX;
+						centerY = tileY * transform_size + transform_size/2 - shiftY;
+						double [][] centersXY = geometryCorrection.getPortsCoordinates(
+								centerX,
+								centerY,
+								disparity);
+						if ((globalDebugLevel > -1) && (tileX >= debug_tileX - 2) && (tileX <= debug_tileX + 2) &&
+								(tileY >= debug_tileY - 2) && (tileY <= debug_tileY+2)) {
+							for (int i = 0; i < quad; i++) {
+							System.out.println("clt_aberrations_quad(): color="+chn+", tileX="+tileX+", tileY="+tileY+
+									" centerX="+centerX+" centerY="+centerY+" disparity="+disparity+
+									" centersXY["+i+"][0]="+centersXY[i][0]+" centersXY["+i+"][1]="+centersXY[i][1]);
+							}
+						}
+						
+						for (int i = 0; i < quad; i++) {
+							fract_shiftsXY[i] = extract_correct_tile( // return a pair of resudual offsets
+									imade_data[i],
+									width,       // image width
+									clt_kernels[i], // [color][tileY][tileX][band][pixel]
+									clt_data[i][chn][tileY][tileX], //double  [][]        clt_tile,    // should be double [4][];
+									kernel_step,
+									transform_size,
+									dtt, 
+									chn,                              
+									centersXY[i][0], // centerX, // center of aberration-corrected (common model) tile, X
+									centersXY[i][1], // centerY, //
+									(globalDebugLevel > 0) && (tileX == debug_tileX) && (tileY == debug_tileY) && (chn == 2), // external tile compare
+									no_deconvolution,
+									transpose);
+						}
+						if ((globalDebugLevel > 0) && (debug_tileX == tileX) && (debug_tileY == tileY)  && (chn == 2)) {
+							showDoubleFloatArrays sdfa_instance = new showDoubleFloatArrays(); // just for debugging?
+							String [] titles = {"CC0","SC0","CS0","SS0","CC1","SC1","CS1","SS1","CC2","SC2","CS2","SS2","CC3","SC3","CS3","SS3"};
+							double [][] dbg_tile = new double [16][];
+							for (int i = 0; i < 16; i++) dbg_tile[i]=clt_data[i>>2][chn][tileY][tileX][i & 3];   
+							sdfa_instance.showArrays(dbg_tile,  transform_size, transform_size, true, "pre-shifted_x"+tileX+"_y"+tileY, titles);
+						}
+						
+						if ((globalDebugLevel > -1) && (tileX >= debug_tileX - 2) && (tileX <= debug_tileX + 2) &&
+								(tileY >= debug_tileY - 2) && (tileY <= debug_tileY+2)) {
+							for (int i = 0; i < quad; i++) {
+							System.out.println("clt_aberrations_quad(): color="+chn+", tileX="+tileX+", tileY="+tileY+
+									" fract_shiftsXY["+i+"][0]="+fract_shiftsXY[i][0]+" fract_shiftsXY["+i+"][1]="+fract_shiftsXY[i][1]);
+							}
+						}
+						
+						if (!no_fract_shift) {
+							// apply residual shift
+							for (int i = 0; i < quad; i++) {
+								fract_shift(    // fractional shift in transform domain. Currently uses sin/cos - change to tables with 2? rotations
+										clt_data[i][chn][tileY][tileX], // double  [][]  clt_tile,
+										transform_size,
+										fract_shiftsXY[i][0],            // double        shiftX,
+										fract_shiftsXY[i][1],            // double        shiftY,
+										//									(globalDebugLevel > 0) && (tileX == debug_tileX) && (tileY == debug_tileY)); // external tile compare
+										((globalDebugLevel > 0) && (chn==0) && (tileX >= debug_tileX - 2) && (tileX <= debug_tileX + 2) &&
+												(tileY >= debug_tileY - 2) && (tileY <= debug_tileY+2)));									
+							}
+							if ((globalDebugLevel > 0) && (debug_tileX == tileX) && (debug_tileY == tileY)) {
+								showDoubleFloatArrays sdfa_instance = new showDoubleFloatArrays(); // just for debugging?
+								String [] titles = {"CC0","SC0","CS0","SS0","CC1","SC1","CS1","SS1","CC2","SC2","CS2","SS2","CC3","SC3","CS3","SS3"};
+								double [][] dbg_tile = new double [16][];
+								for (int i = 0; i < 16; i++) dbg_tile[i]=clt_data[i>>2][chn][tileY][tileX][i & 3];   
+								sdfa_instance.showArrays(dbg_tile,  transform_size, transform_size, true, "shifted_x"+tileX+"_y"+tileY, titles);
+							}
+						}
+					}
+				}
+			};
+		}		      
+		startAndJoin(threads);
+		return clt_data;
+	}
+
+	
+	
 	
 	
 	

src/main/java/PixelMapping.java

@@ -15584,7 +15584,13 @@ public class PixelMapping {
 	    public int    subcamera= -1;
 	    public int    subchannel=-1;
 	    public int    sensor_port=-1;
-	    // TODO: add serial# (and temperature?)
+
+// Following 4 parameters are not included in the vector, names, descriptions
+	    public boolean cartesian = false; // cartesian coordinates mode (false - cylindrical)
+    	public double right;   // distance to the right (radius*sin(azimuth))
+    	public double forward; // distance forward (radius*cos(azimuth))
+    	public double heading;   // absolute heading in degrees (used in cartesian mode)
+	    
     	public double azimuth; // azimuth of the lens entrance pupil center, degrees, clockwise looking from top
     	public double radius;  // mm, distance from the rotation axis
     	public double height;       // mm, up - from the origin point
@@ -16266,8 +16272,33 @@ public class PixelMapping {
 				if (imp.getProperty("r_od_"+i+"_o")  !=null) this.r_od[i][0]= Double.parseDouble((String) imp.getProperty("r_od_"+i+"_o"));
 				if (imp.getProperty("r_od_"+i+"_d")  !=null) this.r_od[i][1]= Double.parseDouble((String) imp.getProperty("r_od_"+i+"_d"));
 			}
-
+			
+        	if (imp.getProperty("forward")  !=null) this.forward=  Double.parseDouble((String) imp.getProperty("forward"));
+        	if (imp.getProperty("right")    !=null) this.right=    Double.parseDouble((String) imp.getProperty("right"));
+        	if (imp.getProperty("aheading") !=null) this.heading=  Double.parseDouble((String) imp.getProperty("aheading"));
+			
+        	if (imp.getProperty("cartesian") !=null) {
+        		this.cartesian= Boolean.parseBoolean((String) imp.getProperty("cartesian"));
+        		updateCartesian(); // recalculate other parameters (they may or may nort be provided
+        	} else {
+        		this.cartesian = false;
+        	}
 	    }
+	    
+    	public void updateCartesian(){ // set alternative parameters
+    		if (cartesian) {
+    			this.azimuth = Math.atan2(this.right, this.forward)*180.0/Math.PI;
+    			this.radius = Math.sqrt(this.forward*this.forward + this.right*this.right);
+    			this.phi = this.heading - this.azimuth; 
+    		} else {
+    			this.forward = this.radius * Math.cos(Math.PI*this.azimuth/180.0); 
+    			this.right =   this.radius * Math.sin(Math.PI*this.azimuth/180.0); 
+    			this.heading = this.phi +  this.azimuth;
+    		}
+    	}
+	    
+	    
+	    
 	    /**
 	     * Calculate rotation matrix that converts sensor coordinates to world coordinates 
 	     * @return 3x3 rotation matrix