Skip to content

I
imagej-elphel
Commit 91616de2 authored by Andrey Filippov's avatar Andrey Filippov
Browse files
Options

trying with actual kernels

parent be7e6c62
Show whitespace changes
Inline Side-by-side
Showing with 508 additions and 217 deletions
src/main/java/EyesisCorrectionParameters.java
View file @ 91616de2
...@@ -1666,7 +1666,15 @@ public class EyesisCorrectionParameters { ...@@ -1666,7 +1666,15 @@ public class EyesisCorrectionParameters {
public double dbg_sigma =2.0; public double dbg_sigma =2.0;
public String dbg_mask = ".........:::::::::.........:::::::::......*..:::::*:::.........:::::::::........."; public String dbg_mask = ".........:::::::::.........:::::::::......*..:::::*:::.........:::::::::.........";
public int dbg_mode = 1; // 0 - old LMA, 1 - new LMA public int dbg_mode = 1; // 0 - old LMA, 1 - new LMA
public int dbg_window_mode = 2; // 0 - none, 1 - square, 2 - sin public int dbg_window_mode = 2; // 0 - none, 1 - square, 2 - sin 3 - sin^2
public boolean centerWindowToTarget = true;
// parameters to extract a kernel from the kernel image file
public int color_channel = 2; // green (<0 - use simulated kernel, also will use simulated if kernels are not set)
public int decimation = 2; // decimate original kernel this much in each direction
public double decimateSigma = 0.4; // what is the optimal value for each decimation?
public int tileX = 82; // number of kernel tile (0..163)
public int tileY = 62; // number of kernel tile (0..122)
public DCTParameters( public DCTParameters(
int dct_size, int dct_size,
...@@ -1697,7 +1705,6 @@ public class EyesisCorrectionParameters { ...@@ -1697,7 +1705,6 @@ public class EyesisCorrectionParameters {
properties.setProperty(prefix+"asym_tax_free", this.asym_tax_free+""); properties.setProperty(prefix+"asym_tax_free", this.asym_tax_free+"");
properties.setProperty(prefix+"seed_size", this.seed_size+""); properties.setProperty(prefix+"seed_size", this.seed_size+"");
properties.setProperty(prefix+"asym_random", this.asym_random+""); properties.setProperty(prefix+"asym_random", this.asym_random+"");
properties.setProperty(prefix+"LMA_steps", this.LMA_steps+""); properties.setProperty(prefix+"LMA_steps", this.LMA_steps+"");
properties.setProperty(prefix+"dbg_x", this.dbg_x+""); properties.setProperty(prefix+"dbg_x", this.dbg_x+"");
properties.setProperty(prefix+"dbg_y", this.dbg_y+""); properties.setProperty(prefix+"dbg_y", this.dbg_y+"");
...@@ -1707,6 +1714,12 @@ public class EyesisCorrectionParameters { ...@@ -1707,6 +1714,12 @@ public class EyesisCorrectionParameters {
properties.setProperty(prefix+"dbg_mask", this.dbg_mask+""); properties.setProperty(prefix+"dbg_mask", this.dbg_mask+"");
properties.setProperty(prefix+"dbg_mode", this.dbg_mode+""); properties.setProperty(prefix+"dbg_mode", this.dbg_mode+"");
properties.setProperty(prefix+"dbg_window_mode", this.dbg_window_mode+""); properties.setProperty(prefix+"dbg_window_mode", this.dbg_window_mode+"");
properties.setProperty(prefix+"centerWindowToTarget", this.centerWindowToTarget+"");
properties.setProperty(prefix+"color_channel", this.color_channel+"");
properties.setProperty(prefix+"decimation", this.dbg_window_mode+"");
properties.setProperty(prefix+"decimateSigma", this.decimateSigma+"");
properties.setProperty(prefix+"tileX", this.tileX+"");
properties.setProperty(prefix+"tileY", this.tileY+"");
} }
public void getProperties(String prefix,Properties properties){ public void getProperties(String prefix,Properties properties){
...@@ -1728,32 +1741,46 @@ public class EyesisCorrectionParameters { ...@@ -1728,32 +1741,46 @@ public class EyesisCorrectionParameters {
if (properties.getProperty(prefix+"dbg_sigma")!=null) this.dbg_sigma=Double.parseDouble(properties.getProperty(prefix+"dbg_sigma")); if (properties.getProperty(prefix+"dbg_sigma")!=null) this.dbg_sigma=Double.parseDouble(properties.getProperty(prefix+"dbg_sigma"));
if (properties.getProperty(prefix+"dbg_mask")!=null) this.dbg_mask=properties.getProperty(prefix+"dbg_mask"); if (properties.getProperty(prefix+"dbg_mask")!=null) this.dbg_mask=properties.getProperty(prefix+"dbg_mask");
if (properties.getProperty(prefix+"dbg_mode")!=null) this.dbg_mode=Integer.parseInt(properties.getProperty(prefix+"dbg_mode")); if (properties.getProperty(prefix+"dbg_mode")!=null) this.dbg_mode=Integer.parseInt(properties.getProperty(prefix+"dbg_mode"));
if (properties.getProperty(prefix+"tileY")!=null) this.tileY=Integer.parseInt(properties.getProperty(prefix+"tileY"));
if (properties.getProperty(prefix+"centerWindowToTarget")!=null) this.centerWindowToTarget=Boolean.parseBoolean(properties.getProperty(prefix+"centerWindowToTarget"));
if (properties.getProperty(prefix+"color_channel")!=null) this.color_channel=Integer.parseInt(properties.getProperty(prefix+"color_channel"));
if (properties.getProperty(prefix+"decimation")!=null) this.decimation=Integer.parseInt(properties.getProperty(prefix+"decimation"));
if (properties.getProperty(prefix+"decimateSigma")!=null) this.decimateSigma=Double.parseDouble(properties.getProperty(prefix+"decimateSigma"));
if (properties.getProperty(prefix+"tileX")!=null) this.tileX=Integer.parseInt(properties.getProperty(prefix+"tileX"));
if (properties.getProperty(prefix+"dbg_window_mode")!=null) this.dbg_window_mode=Integer.parseInt(properties.getProperty(prefix+"dbg_window_mode")); if (properties.getProperty(prefix+"dbg_window_mode")!=null) this.dbg_window_mode=Integer.parseInt(properties.getProperty(prefix+"dbg_window_mode"));
} }
public boolean showDialog() { public boolean showDialog() {
GenericDialog gd = new GenericDialog("Set DCT parameters"); GenericDialog gd = new GenericDialog("Set DCT parameters");
gd.addNumericField("DCT size", this.dct_size, 0); //32 gd.addNumericField("DCT size", this.dct_size, 0);
gd.addNumericField("Size of asymmetrical (non-DCT) kernel", this.asym_size, 0); //6 gd.addNumericField("Size of asymmetrical (non-DCT) kernel", this.asym_size, 0);
gd.addNumericField("Maximal number of non-zero pixels in direct convolution kernel", this.asym_pixels, 0); //6 gd.addNumericField("Maximal number of non-zero pixels in direct convolution kernel", this.asym_pixels, 0);
gd.addNumericField("How far to try a new asym kernel pixel from existing ones", this.asym_distance, 0); //6 gd.addNumericField("How far to try a new asym kernel pixel from existing ones", this.asym_distance, 0);
gd.addNumericField("MDCT window type (0,1,2)", this.dct_window, 0); //0..2 gd.addNumericField("MDCT window type (0,1,2)", this.dct_window, 0);
gd.addNumericField("LMA_steps", this.LMA_steps, 0); //0..2 gd.addNumericField("LMA_steps", this.LMA_steps, 0);
gd.addNumericField("Compactness (punish off-center asym_kernel pixels (proportional to r^2)", this.compactness, 2); //0..2 gd.addNumericField("Compactness (punish off-center asym_kernel pixels (proportional to r^2)", this.compactness,2);
gd.addNumericField("Factorization target precision (stop if achieved)", this.fact_precision, 4); //0..2 gd.addNumericField("Factorization target precision (stop if achieved)", this.fact_precision, 4);
gd.addNumericField("Do not punish pixels in the square around center", this.asym_tax_free, 0); //0..2 gd.addNumericField("Do not punish pixels in the square around center", this.asym_tax_free, 0);
gd.addNumericField("Start asym_kernel with this number of pixels (0 - single, 4n+0 (X between cells), 4*n+1 - x around center cell", this.seed_size, 0); //0..2 gd.addNumericField("Start asym_kernel with this number of pixels (0 - single, 4n+0 (X between cells), 4*n+1 - x around center cell", this.seed_size, 0); //0..2
gd.addNumericField("Initialize asym_kernel with random numbers (amplitude)", this.asym_random, 2); //0..2 gd.addNumericField("Initialize asym_kernel with random numbers (amplitude)", this.asym_random, 2);
gd.addNumericField("dbg_x", this.dbg_x, 2); //0..2 gd.addNumericField("dbg_x", this.dbg_x, 2);
gd.addNumericField("dbg_y", this.dbg_y, 2); //0..2 gd.addNumericField("dbg_y", this.dbg_y, 2);
gd.addNumericField("dbg_x1", this.dbg_x1, 2); //0..2 gd.addNumericField("dbg_x1", this.dbg_x1, 2);
gd.addNumericField("dbg_y1", this.dbg_y1, 2); //0..2 gd.addNumericField("dbg_y1", this.dbg_y1, 2);
gd.addNumericField("dbg_sigma", this.dbg_sigma, 3); //0..2 gd.addNumericField("dbg_sigma", this.dbg_sigma, 3);
gd.addStringField ("Debug mask (anything but * is false)", this.dbg_mask,100); gd.addStringField ("Debug mask (anything but * is false)", this.dbg_mask, 100);
gd.addNumericField("LMA implementation: 0 - old, 1 - new", this.dbg_mode, 0); //32 gd.addNumericField("LMA implementation: 0 - old, 1 - new", this.dbg_mode, 0);
gd.addNumericField("Convolution window: 0 - none, 1 - square, 2 - sin, 3 - sin^2", this.dbg_window_mode, 0); //32 gd.addNumericField("Convolution window: 0 - none, 1 - square, 2 - sin, 3 - sin^2", this.dbg_window_mode, 0);
gd.addCheckbox ("Center convolution window around target kernel center", this.centerWindowToTarget);
// gd.addNumericField("Debug Level:", MASTER_DEBUG_LEVEL, 0); gd.addNumericField("Color channel to extract kernel (<0 - use synthetic)", this.color_channel, 0);
gd.addNumericField("Convolution kernel decimation (original is normally 2x)", this.decimation, 0);
gd.addNumericField("Smooth convolution kernel before decimation", this.decimateSigma, 3);
gd.addNumericField("Tile X to extract (0..163)", this.tileX, 0);
gd.addNumericField("Tile Y to extract (0..122)", this.tileY, 0);
gd.showDialog(); gd.showDialog();
if (gd.wasCanceled()) return false; if (gd.wasCanceled()) return false;
this.dct_size= (int) gd.getNextNumber(); this.dct_size= (int) gd.getNextNumber();
this.asym_size= (int) gd.getNextNumber(); this.asym_size= (int) gd.getNextNumber();
...@@ -1774,6 +1801,12 @@ public class EyesisCorrectionParameters { ...@@ -1774,6 +1801,12 @@ public class EyesisCorrectionParameters {
this.dbg_mask= gd.getNextString(); this.dbg_mask= gd.getNextString();
this.dbg_mode= (int) gd.getNextNumber(); this.dbg_mode= (int) gd.getNextNumber();
this.dbg_window_mode= (int) gd.getNextNumber(); this.dbg_window_mode= (int) gd.getNextNumber();
this.centerWindowToTarget= gd.getNextBoolean();
this.color_channel= (int) gd.getNextNumber();
this.decimation= (int) gd.getNextNumber();
this.decimateSigma= gd.getNextNumber();
this.tileX= (int) gd.getNextNumber();
this.tileY= (int) gd.getNextNumber();
// MASTER_DEBUG_LEVEL= (int) gd.getNextNumber(); // MASTER_DEBUG_LEVEL= (int) gd.getNextNumber();
return true; return true;
......
src/main/java/EyesisDCT.java
View file @ 91616de2
...@@ -32,7 +32,8 @@ public class EyesisDCT { ...@@ -32,7 +32,8 @@ public class EyesisDCT {
public EyesisCorrections eyesisCorrections = null; public EyesisCorrections eyesisCorrections = null;
public EyesisCorrectionParameters.CorrectionParameters correctionsParameters=null; public EyesisCorrectionParameters.CorrectionParameters correctionsParameters=null;
public EyesisCorrectionParameters.DCTParameters dctParameters = null; public EyesisCorrectionParameters.DCTParameters dctParameters = null;
public DCTKernel [] kernels = null; public DCTKernels [] kernels = null;
public ImagePlus eyesisKernelImage = null;
public EyesisDCT( public EyesisDCT(
EyesisCorrections eyesisCorrections, EyesisCorrections eyesisCorrections,
EyesisCorrectionParameters.CorrectionParameters correctionsParameters, EyesisCorrectionParameters.CorrectionParameters correctionsParameters,
...@@ -42,13 +43,35 @@ public class EyesisDCT { ...@@ -42,13 +43,35 @@ public class EyesisDCT {
this.correctionsParameters = correctionsParameters; this.correctionsParameters = correctionsParameters;
this.dctParameters= dctParameters; this.dctParameters= dctParameters;
} }
public class DCTKernel{ public class DCTKernels{
// -- old --
public int size = 32; // kernel (DCT) size public int size = 32; // kernel (DCT) size
public int img_step = 32; // pixel step in the image for each kernel public int img_step = 32; // pixel step in the image for each kernel
public double [][][] offsets = null; // per color, per kernel,per coord public double [][][] offsets = null; // per color, per kernel,per coord
public double [][] kern = null; // kernel image in linescan order public double [][] kern = null; // kernel image in linescan order
// -- new --
public int numHor = 164; // number of kernel tiles in a row
public int dct_size = 8; // DCT-II size, sym. kernel square side is 2*dct_size-1
public int asym_size = 15; // asymmetrical convolution limits, odd
public int asym_nonzdero = 4; // maximal number of non-zero elements in the asymmetrical kernels
public double [][] sym_kernels = null; // per-color channel, DCT kernels in linescan order
public double [][] asym_kernels = null; // per-color channel, asymmetrical kernels (all but asym_nonzdero elements are strictly 0)
} }
public void setKernelImageFile(ImagePlus img_kernels){
eyesisKernelImage = img_kernels;
}
public boolean kernelImageSet(){
return eyesisKernelImage != null;
}
public boolean createDCTKernels( public boolean createDCTKernels(
EyesisCorrectionParameters.DCTParameters dct_parameters,
int srcKernelSize,
int threadsMax, // maximal number of threads to launch
boolean updateStatus,
int debugLevel
){ ){
String [] sharpKernelPaths= correctionsParameters.selectKernelChannelFiles( String [] sharpKernelPaths= correctionsParameters.selectKernelChannelFiles(
0, // 0 - sharp, 1 - smooth 0, // 0 - sharp, 1 - smooth
...@@ -59,11 +82,13 @@ public class EyesisDCT { ...@@ -59,11 +82,13 @@ public class EyesisDCT {
System.out.println(i+":"+sharpKernelPaths[i]); System.out.println(i+":"+sharpKernelPaths[i]);
} }
if (kernels == null){ if (kernels == null){
kernels = new DCTKernel[eyesisCorrections.usedChannels.length]; kernels = new DCTKernels[eyesisCorrections.usedChannels.length];
for (int chn=0;chn<eyesisCorrections.usedChannels.length;chn++){ for (int chn=0;chn<eyesisCorrections.usedChannels.length;chn++){
kernels[chn] = null; kernels[chn] = null;
} }
} }
showDoubleFloatArrays sdfa_instance = new showDoubleFloatArrays(); // just for debugging?
for (int chn=0;chn<eyesisCorrections.usedChannels.length;chn++){ for (int chn=0;chn<eyesisCorrections.usedChannels.length;chn++){
if (eyesisCorrections.usedChannels[chn] && (sharpKernelPaths[chn]!=null) && (kernels[chn]==null)){ if (eyesisCorrections.usedChannels[chn] && (sharpKernelPaths[chn]!=null) && (kernels[chn]==null)){
ImagePlus imp_kernel_sharp=new ImagePlus(sharpKernelPaths[chn]); ImagePlus imp_kernel_sharp=new ImagePlus(sharpKernelPaths[chn]);
...@@ -73,24 +98,33 @@ public class EyesisDCT { ...@@ -73,24 +98,33 @@ public class EyesisDCT {
continue; continue;
} }
ImageStack kernel_sharp_stack= imp_kernel_sharp.getStack(); ImageStack kernel_sharp_stack= imp_kernel_sharp.getStack();
System.out.println("debugLevel = "+debugLevel+" kernel_sharp_stack.getWidth() = "+kernel_sharp_stack.getWidth()+
" kernel_sharp_stack.getHeight() = "+kernel_sharp_stack.getHeight());
DCTKernels kernels = calculateDCTKernel (
kernel_sharp_stack, // final ImageStack kernelStack, // first stack with 3 colors/slices convolution kernels
srcKernelSize, // final int kernelSize, // 64
dct_parameters, // final double blurSigma,
threadsMax, // maximal number of threads to launch
updateStatus,
debugLevel); // update status info
int sym_width = kernels.numHor * kernels.dct_size;
int sym_height = kernels.sym_kernels[0].length /sym_width;
sdfa_instance.showArrays(kernels.sym_kernels, sym_width, sym_height, true, imp_kernel_sharp.getTitle()+"-sym");
int asym_width = kernels.numHor * kernels.asym_size;
int asym_height = kernels.asym_kernels[0].length /asym_width;
sdfa_instance.showArrays(kernels.asym_kernels, asym_width, asym_height, true, imp_kernel_sharp.getTitle()+"-asym");
} }
} }
return true; return true;
} }
public DCTKernel calculateDCTKernel ( public DCTKernels calculateDCTKernel (
final ImageStack kernelStack, // first stack with 3 colors/slices convolution kernels final ImageStack kernelStack, // first stack with 3 colors/slices convolution kernels
final int kernelSize, // 64 final int kernelSize, // 64
final double blurSigma, final EyesisCorrectionParameters.DCTParameters dct_parameters,
final int scaleDown, // kernels are saved with higher resolution - scale them down by this value (2)
final int dctSize,
final int threadsMax, // maximal number of threads to launch final int threadsMax, // maximal number of threads to launch
final boolean updateStatus, final boolean updateStatus,
final int globalDebugLevel) // update status info final int globalDebugLevel) // update status info
...@@ -100,29 +134,38 @@ public class EyesisDCT { ...@@ -100,29 +134,38 @@ public class EyesisDCT {
final int kernelNumHor=kernelWidth/kernelSize; final int kernelNumHor=kernelWidth/kernelSize;
final int kernelNumVert=kernelStack.getHeight()/kernelSize; final int kernelNumVert=kernelStack.getHeight()/kernelSize;
final int nChn=kernelStack.getSize(); final int nChn=kernelStack.getSize();
final int length=kernelNumHor*kernelNumVert*dctSize*dctSize;// size of kernel data // final int length=kernelNumHor*kernelNumVert* dct_parameters.dct_size * dct_parameters.dct_size;// size of kernel data
final DCTKernel dct_kernel = new DCTKernel(); final DCTKernels dct_kernel = new DCTKernels();
dct_kernel.size = dctSize; dct_kernel.size = dct_parameters.dct_size;
dct_kernel.img_step = dctSize; // configure? dct_kernel.img_step = kernelSize/2/dct_parameters.decimation ; // May be wrong
dct_kernel.offsets = new double[nChn][kernelNumHor*kernelNumVert][2]; dct_kernel.sym_kernels = new double [nChn][kernelNumHor*kernelNumVert*dct_parameters.dct_size * dct_parameters.dct_size];
dct_kernel.kern = new double [nChn][kernelNumHor*kernelNumVert*dctSize*dctSize]; dct_kernel.asym_kernels = new double [nChn][kernelNumHor*kernelNumVert*dct_parameters.asym_size * dct_parameters.asym_size];
// currently each 64x64 kernel corresponds to 16x16 original pixels tile, 2 tiles margin each side // currently each 64x64 kernel corresponds to 16x16 original pixels tile, 2 tiles margin each side
final Thread[] threads = ImageDtt.newThreadArray(threadsMax); final Thread[] threads = ImageDtt.newThreadArray(threadsMax);
final AtomicInteger ai = new AtomicInteger(0); final AtomicInteger ai = new AtomicInteger(0);
final int numberOfKernels= kernelNumHor*kernelNumVert*nChn; final int numberOfKernels= kernelNumHor*kernelNumVert*nChn;
final int numberOfKernelsInChn=kernelNumHor*kernelNumVert; final int numberOfKernelsInChn=kernelNumHor*kernelNumVert;
final long startTime = System.nanoTime(); final long startTime = System.nanoTime();
System.out.println("calculateDCTKernel():numberOfKernels="+numberOfKernels);
for (int ithread = 0; ithread < threads.length; ithread++) { for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() { threads[ithread] = new Thread() {
public void run() { public void run() {
DoubleGaussianBlur gb=null; DoubleGaussianBlur gb=null;
if (blurSigma>0) gb=new DoubleGaussianBlur(); if (dct_parameters.decimateSigma > 0) gb=new DoubleGaussianBlur();
float [] kernelPixels= null; // will be initialized at first use float [] kernelPixels= null; // will be initialized at first use NOT yet?
double [] kernel= new double[kernelSize*kernelSize]; double [] kernel= new double[kernelSize*kernelSize];
int targetSize = dct_parameters.asym_size + 2 * dct_parameters.dct_size - 2;
double [] target_kernel = new double [targetSize * targetSize];
FactorConvKernel factorConvKernel = new FactorConvKernel();
factorConvKernel.setDebugLevel (0); // globalDebugLevel);
factorConvKernel.setTargetWindowMode (dct_parameters.dbg_window_mode, dct_parameters.centerWindowToTarget);
factorConvKernel.numIterations = dct_parameters.LMA_steps;
factorConvKernel.setAsymCompactness (dct_parameters.compactness, dct_parameters.asym_tax_free);
int chn,tileY,tileX; int chn,tileY,tileX;
int chn0=-1; // int chn0=-1;
int i; // int i;
double sum; // double sum;
for (int nTile = ai.getAndIncrement(); nTile < numberOfKernels; nTile = ai.getAndIncrement()) { for (int nTile = ai.getAndIncrement(); nTile < numberOfKernels; nTile = ai.getAndIncrement()) {
chn=nTile/numberOfKernelsInChn; chn=nTile/numberOfKernelsInChn;
tileY =(nTile % numberOfKernelsInChn)/kernelNumHor; tileY =(nTile % numberOfKernelsInChn)/kernelNumHor;
...@@ -131,32 +174,94 @@ public class EyesisDCT { ...@@ -131,32 +174,94 @@ public class EyesisDCT {
if (updateStatus) IJ.showStatus("Processing kernels, channel "+(chn+1)+" of "+nChn+", row "+(tileY+1)+" of "+kernelNumVert); if (updateStatus) IJ.showStatus("Processing kernels, channel "+(chn+1)+" of "+nChn+", row "+(tileY+1)+" of "+kernelNumVert);
if (globalDebugLevel>2) System.out.println("Processing kernels, channel "+(chn+1)+" of "+nChn+", row "+(tileY+1)+" of "+kernelNumVert+" : "+IJ.d2s(0.000000001*(System.nanoTime()-startTime),3)); if (globalDebugLevel>2) System.out.println("Processing kernels, channel "+(chn+1)+" of "+nChn+", row "+(tileY+1)+" of "+kernelNumVert+" : "+IJ.d2s(0.000000001*(System.nanoTime()-startTime),3));
} }
kernelPixels=(float[]) kernelStack.getPixels(chn+1);
/* read convolution kernel */ /* read convolution kernel */
extractOneKernel(kernelPixels, // array of combined square kernels, each extractOneKernel(
kernelPixels, // array of combined square kernels, each
kernel, // will be filled, should have correct size before call kernel, // will be filled, should have correct size before call
kernelNumHor, // number of kernels in a row kernelNumHor, // number of kernels in a row
tileX, // horizontal number of kernel to extract tileX, // horizontal number of kernel to extract
tileY); // vertical number of kernel to extract tileY); // vertical number of kernel to extract
if (blurSigma>0) gb.blurDouble(kernel, kernelSize, kernelSize, blurSigma, blurSigma, 0.01);
/* Calculate sum of squared kernel1 elements */
sum=0.0;
for (i=0; i<kernel.length;i++) sum+=kernel[i]*kernel[i];
// outPixles[chn][tileY*kernelNumHor+tileX]= (float) (Math.sqrt(sum)); reformatKernel(
// System.out.println("Processing kernels, channel "+(chn+1)+" of "+nChn+", row "+(tileY+1)+" of "+kernelNumVert+" sum="+sum); kernel, // will be blured in-place
target_kernel, // expand/crop, blur/decimate result
kernelSize,
targetSize,
dct_parameters.decimation,
dct_parameters.decimateSigma,
gb);
// int numAsym =
factorConvKernel.calcKernels(
target_kernel,
dct_parameters.asym_size,
dct_parameters.dct_size,
dct_parameters.fact_precision,
dct_parameters.asym_pixels, // maximal number of non-zero pixels in asymmmetrical kernel
dct_parameters.asym_distance, // how far to seed a new pixel
dct_parameters.seed_size);
double [] sym_kernel = factorConvKernel.getSymKernel();
double [] asym_kernel = factorConvKernel.getAsymKernel();
int sym_kernel_inc_index = kernelNumHor * dct_parameters.dct_size;
int sym_kernel_start_index = (sym_kernel_inc_index * tileY + tileX) * dct_parameters.dct_size;
// int indx = 0;
for (int i = 0; i<dct_parameters.dct_size;i++){
System.arraycopy(
sym_kernel,
i * dct_parameters.dct_size,
dct_kernel.sym_kernels[chn],
sym_kernel_start_index + i * sym_kernel_inc_index,
dct_parameters.dct_size);
/*
int dst_start = sym_kernel_start_index + i * sym_kernel_inc_index;
for (int j = 0; j < dct_parameters.dct_size; j++){
dct_kernel.sym_kernels[chn][dst_start++] = sym_kernel[indx++];
}
*/
}
int asym_kernel_inc_index = kernelNumHor * dct_parameters.asym_size;
int asym_kernel_start_index = (asym_kernel_inc_index * tileY + tileX)* dct_parameters.asym_size;
for (int i = 0; i<dct_parameters.asym_size;i++){
System.arraycopy(
asym_kernel,
i * dct_parameters.asym_size,
dct_kernel.asym_kernels[chn],
asym_kernel_start_index + i * asym_kernel_inc_index,
dct_parameters.asym_size);
}
} }
} }
}; };
} }
ImageDtt.startAndJoin(threads); ImageDtt.startAndJoin(threads);
if (globalDebugLevel > 1) System.out.println("Threads done at "+IJ.d2s(0.000000001*(System.nanoTime()-startTime),3)); if (globalDebugLevel > 1) System.out.println("Threads done at "+IJ.d2s(0.000000001*(System.nanoTime()-startTime),3));
System.out.println("1.Threads done at "+IJ.d2s(0.000000001*(System.nanoTime()-startTime),3));
/* prepare result stack to return */ /* prepare result stack to return */
ImageStack outStack=new ImageStack(kernelNumHor,kernelNumVert); // ImageStack outStack=new ImageStack(kernelNumHor,kernelNumVert);
return dct_kernel; return dct_kernel;
} }
/*
* final int kernelNumHor=kernelWidth/kernelSize;
final int kernelNumVert=kernelStack.getHeight()/kernelSize;
final int nChn=kernelStack.getSize();
dct_kernel.sym_kernels = new double [nChn][kernelNumHor*kernelNumVert*dct_parameters.dct_size * dct_parameters.dct_size];
dct_kernel.asym_kernels = new double [nChn][kernelNumHor*kernelNumVert*dct_parameters.asym_size * dct_parameters.asym_size];
*
System.arraycopy(currentfX, 0, convolved, 0, convolved.length);
DCT_PARAMETERS.fact_precision,
DCT_PARAMETERS.asym_pixels,
DCT_PARAMETERS.asym_distance,
DCT_PARAMETERS.seed_size);
*/
//processChannelImage //processChannelImage
//convolveStackWithKernelStack //convolveStackWithKernelStack
// Remove later, copied here as a reference
public ImageStack calculateKernelsNoiseGains ( public ImageStack calculateKernelsNoiseGains (
final ImageStack kernelStack1, // first stack with 3 colors/slices convolution kernels final ImageStack kernelStack1, // first stack with 3 colors/slices convolution kernels
final ImageStack kernelStack2, // second stack with 3 colors/slices convolution kernels (or null) final ImageStack kernelStack2, // second stack with 3 colors/slices convolution kernels (or null)
...@@ -243,7 +348,48 @@ public class EyesisDCT { ...@@ -243,7 +348,48 @@ public class EyesisDCT {
return outStack; return outStack;
} }
void extractOneKernel(float [] pixels, // array of combined square kernels, each // mostly for testing
//eyesisKernelImage
public double [] extractOneKernelFromStack(
final int kernelSize, // 64
final int chn,
final int xTile, // horizontal number of kernel to extract
final int yTile) // vertical number of kernel to extract
{
if (eyesisKernelImage == null) return null;
final ImageStack kernelStack = eyesisKernelImage.getStack();
return extractOneKernelFromStack(
kernelStack, // first stack with 3 colors/slices convolution kernels
kernelSize, // 64
chn,
xTile, // horizontal number of kernel to extract
yTile); // vertical number of kernel to extract
}
public double [] extractOneKernelFromStack(
final ImageStack kernelStack, // first stack with 3 colors/slices convolution kernels
final int kernelSize, // 64
final int chn,
final int xTile, // horizontal number of kernel to extract
final int yTile) // vertical number of kernel to extract
{
final int kernelWidth=kernelStack.getWidth();
final int kernelNumHor=kernelWidth/kernelSize;
// final int kernelNumVert=kernelStack.getHeight()/kernelSize;
// final int nChn=kernelStack.getSize();
double [] kernel = new double [kernelSize*kernelSize];
extractOneKernel((float[]) kernelStack.getPixels(chn+1), // array of combined square kernels, each
kernel, // will be filled, should have correct size before call
kernelNumHor, // number of kernels in a row
xTile, // horizontal number of kernel to extract
yTile);
return kernel;
}
//imp2.getStack()
// to be used in threaded method
private void extractOneKernel(float [] pixels, // array of combined square kernels, each
double [] kernel, // will be filled, should have correct size before call double [] kernel, // will be filled, should have correct size before call
int numHor, // number of kernels in a row int numHor, // number of kernels in a row
int xTile, // horizontal number of kernel to extract int xTile, // horizontal number of kernel to extract
...@@ -262,8 +408,55 @@ public class EyesisDCT { ...@@ -262,8 +408,55 @@ public class EyesisDCT {
int base=(yTile*pixelsWidth+xTile)*size; int base=(yTile*pixelsWidth+xTile)*size;
for (i=0;i<size;i++) for (j=0;j<size;j++) kernel [i*size+j]=pixels[base+i*pixelsWidth+j]; for (i=0;i<size;i++) for (j=0;j<size;j++) kernel [i*size+j]=pixels[base+i*pixelsWidth+j];
} }
public double [] reformatKernel(
double [] src_kernel,// will be blured in-place
int src_size, // typical 64
int dst_size, // typical 15
int decimation,// typical 2
double sigma)
{
double [] dst_kernel = new double [dst_size*dst_size];
DoubleGaussianBlur gb = null;
if (sigma > 0) gb = new DoubleGaussianBlur();
reformatKernel(
src_kernel,
dst_kernel,
src_size,
dst_size,
decimation,
sigma,
gb);
return dst_kernel;
}
// to be used in threaded method
private void reformatKernel(
double [] src_kernel, // will be blured in-place
double [] dst_kernel,
int src_size,
int dst_size,
int decimation,
double sigma,
DoubleGaussianBlur gb)
{
if (gb != null) gb.blurDouble(src_kernel, src_size, src_size, sigma, sigma, 0.01);
int src_center = src_size / 2; // 32
int dst_center = dst_size / 2; // 7
for (int i = 0; i< dst_size; i++){
int src_i = (i - dst_center)*decimation + src_center;
if ((src_i >= 0) && (src_i < src_size)) {
for (int j = 0; j< dst_size; j++) {
int src_j = (j - dst_center)*decimation + src_center;
if ((src_j >= 0) && (src_j < src_size)) {
dst_kernel[i*dst_size + j] = src_kernel[src_i*src_size + src_j];
} else {
dst_kernel[i*dst_size + j] = 0;
}
}
} else {
for (int j = 0; j< dst_size; j++) dst_kernel[i*dst_size + j] = 0;
}
}
}
} }
src/main/java/Eyesis_Correction.java
View file @ 91616de2
...@@ -450,13 +450,14 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -450,13 +450,14 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
if (DCT_MODE) { if (DCT_MODE) {
panelDct1 = new Panel(); panelDct1 = new Panel();
panelDct1.setLayout(new GridLayout(1, 0, 5, 5)); // rows, columns, vgap, hgap panelDct1.setLayout(new GridLayout(1, 0, 5, 5)); // rows, columns, vgap, hgap
addButton("DCT test 1",panelDct1,color_process); addButton("DCT test 1", panelDct1, color_process);
addButton("DCT test 2",panelDct1,color_process); addButton("DCT test 2", panelDct1, color_process);
addButton("DCT test 3",panelDct1,color_process); addButton("DCT test 3", panelDct1, color_process);
addButton("DCT test 4",panelDct1,color_process); addButton("DCT test 4", panelDct1, color_process);
addButton("Test Kernel Factorization",panelDct1,color_process); addButton("Test Kernel Factorization", panelDct1, color_process);
addButton("Min Kernel Factorization",panelDct1,color_process); addButton("Min Kernel Factorization", panelDct1, color_process);
addButton("Create DCT kernels",panelDct1,color_process); addButton("Select kernels image", panelDct1, color_configure);
addButton("Create DCT kernels", panelDct1, color_process);
add(panelDct1); add(panelDct1);
} }
pack(); pack();
...@@ -2838,14 +2839,33 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2838,14 +2839,33 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
if (!DCT_PARAMETERS.showDialog()) return; if (!DCT_PARAMETERS.showDialog()) return;
FactorConvKernel factorConvKernel = new FactorConvKernel(DCT_PARAMETERS.dbg_mode == 1); FactorConvKernel factorConvKernel = new FactorConvKernel(DCT_PARAMETERS.dbg_mode == 1);
factorConvKernel.setDebugLevel(DEBUG_LEVEL); factorConvKernel.setDebugLevel(DEBUG_LEVEL);
factorConvKernel.setTargetWindowMode(DCT_PARAMETERS.dbg_window_mode); factorConvKernel.setTargetWindowMode(DCT_PARAMETERS.dbg_window_mode, DCT_PARAMETERS.centerWindowToTarget);
factorConvKernel.numIterations = DCT_PARAMETERS.LMA_steps; factorConvKernel.numIterations = DCT_PARAMETERS.LMA_steps;
factorConvKernel.setAsymCompactness( factorConvKernel.setAsymCompactness(
DCT_PARAMETERS.compactness, DCT_PARAMETERS.compactness,
DCT_PARAMETERS.asym_tax_free); DCT_PARAMETERS.asym_tax_free);
int target_kernel_size = 2*DCT_PARAMETERS.dct_size - 1; int target_kernel_size = 2*DCT_PARAMETERS.dct_size - 1;
int target_expanded_size = 2*DCT_PARAMETERS.dct_size + DCT_PARAMETERS.asym_size -2;
double [] target_expanded = null;
if ((EYESIS_DCT != null) && EYESIS_DCT.kernelImageSet() && (DCT_PARAMETERS.color_channel >= 0)){
System.out.println("Using extracted target kernel");
double [] src_kernel = EYESIS_DCT.extractOneKernelFromStack(
CONVOLVE_FFT_SIZE/2 , // 64
DCT_PARAMETERS.color_channel, // 0..2
DCT_PARAMETERS.tileX, // horizontal number of kernel to extract
DCT_PARAMETERS.tileY); // vertical number of kernel to extract
target_expanded = EYESIS_DCT.reformatKernel(
src_kernel,// will be blured in-place
CONVOLVE_FFT_SIZE/2, // typical 64
target_expanded_size, // typical 15
DCT_PARAMETERS.decimation,// typical 2
DCT_PARAMETERS.decimateSigma);
} else {
System.out.println("Using synthesized target kernel");
double [] target_kernel = new double [target_kernel_size * target_kernel_size]; double [] target_kernel = new double [target_kernel_size * target_kernel_size];
for (int ii=0; ii < target_kernel.length; ii++) target_kernel[ii]=0.0; for (int ii=0; ii < target_kernel.length; ii++) target_kernel[ii]=0.0;
double dist = Math.sqrt((DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)*(DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)+ double dist = Math.sqrt((DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)*(DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)+
(DCT_PARAMETERS.dbg_y1-DCT_PARAMETERS.dbg_y)*(DCT_PARAMETERS.dbg_y1-DCT_PARAMETERS.dbg_y)); (DCT_PARAMETERS.dbg_y1-DCT_PARAMETERS.dbg_y)*(DCT_PARAMETERS.dbg_y1-DCT_PARAMETERS.dbg_y));
...@@ -2866,7 +2886,7 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2866,7 +2886,7 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
DoubleGaussianBlur gb=null; DoubleGaussianBlur gb=null;
if (blurSigma>0) gb=new DoubleGaussianBlur(); if (blurSigma>0) gb=new DoubleGaussianBlur();
if (blurSigma>0) gb.blurDouble(target_kernel, target_kernel_size, target_kernel_size, blurSigma, blurSigma, 0.01); if (blurSigma>0) gb.blurDouble(target_kernel, target_kernel_size, target_kernel_size, blurSigma, blurSigma, 0.01);
// SDFA_INSTANCE.showArrays(target_kernel, target_kernel_size, target_kernel_size, "target_kernel"); // SDFA_INSTANCE.showArrays(target_kernel, target_kernel_size, target_kernel_size, "target_kernel");
if ( if (
(DCT_PARAMETERS.dbg_x == 0.0) && (DCT_PARAMETERS.dbg_x == 0.0) &&
...@@ -2880,10 +2900,8 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2880,10 +2900,8 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
target_kernel[target_kernel_size*ii+jj] = target_kernel[target_kernel_size*ii1+jj1]; target_kernel[target_kernel_size*ii+jj] = target_kernel[target_kernel_size*ii1+jj1];
} }
} }
target_expanded = new double [target_expanded_size * target_expanded_size];
int target_expanded_size = 2*DCT_PARAMETERS.dct_size + DCT_PARAMETERS.asym_size -2;
double [] target_expanded = new double [target_expanded_size * target_expanded_size];
for (int ii=0; ii < target_expanded.length; ii++) target_expanded[ii]=0.0; for (int ii=0; ii < target_expanded.length; ii++) target_expanded[ii]=0.0;
int left_top_margin = ((DCT_PARAMETERS.asym_size-1)/2); int left_top_margin = ((DCT_PARAMETERS.asym_size-1)/2);
for (int ii=0;ii < target_kernel_size; ii++){ for (int ii=0;ii < target_kernel_size; ii++){
...@@ -2892,6 +2910,7 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2892,6 +2910,7 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
target_kernel[ii*target_kernel_size + jj]; target_kernel[ii*target_kernel_size + jj];
} }
} }
}
boolean [] mask = null; boolean [] mask = null;
if (!DCT_PARAMETERS.dbg_mask.equals("")){ if (!DCT_PARAMETERS.dbg_mask.equals("")){
...@@ -2901,7 +2920,7 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2901,7 +2920,7 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
} }
} }
long startTime = System.nanoTime();
boolean result = factorConvKernel.calcKernels( boolean result = factorConvKernel.calcKernels(
target_expanded, target_expanded,
DCT_PARAMETERS.asym_size, DCT_PARAMETERS.asym_size,
...@@ -2911,6 +2930,11 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2911,6 +2930,11 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
DCT_PARAMETERS.seed_size, DCT_PARAMETERS.seed_size,
DCT_PARAMETERS.asym_random); DCT_PARAMETERS.asym_random);
System.out.println("factorConvKernel.calcKernels() returned: >>> "+result+ " <<<"); System.out.println("factorConvKernel.calcKernels() returned: >>> "+result+ " <<<");
System.out.println(
"RMS = "+factorConvKernel.getRMSes()[0]+
", RMSPure = "+factorConvKernel.getRMSes()[1]+
", relRMSPure = "+(factorConvKernel.getRMSes()[1]/factorConvKernel.getTargetRMS())+
", spent "+ IJ.d2s(0.000000001*(System.nanoTime()-startTime),3)+" sec");
double [] sym_kernel = factorConvKernel.getSymKernel(); double [] sym_kernel = factorConvKernel.getSymKernel();
double [] asym_kernel = factorConvKernel.getAsymKernel(); double [] asym_kernel = factorConvKernel.getAsymKernel();
double [] convolved = factorConvKernel.getConvolved(); double [] convolved = factorConvKernel.getConvolved();
...@@ -2923,11 +2947,12 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2923,11 +2947,12 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
weighted_diff100[ii] = diff100[ii]* target_weights[ii]; weighted_diff100[ii] = diff100[ii]* target_weights[ii];
} }
double [][] compare_kernels = {target_expanded, convolved, weighted_diff100,target_weights, diff100}; double [][] compare_kernels = {target_expanded, convolved, weighted_diff100,target_weights, diff100};
if (DEBUG_LEVEL>0) {
System.out.println("DCT_PARAMETERS.dct_size="+DCT_PARAMETERS.dct_size+" DCT_PARAMETERS.asym_size="+DCT_PARAMETERS.asym_size); System.out.println("DCT_PARAMETERS.dct_size="+DCT_PARAMETERS.dct_size+" DCT_PARAMETERS.asym_size="+DCT_PARAMETERS.asym_size);
System.out.println("sym_kernel.length="+ sym_kernel.length); System.out.println("sym_kernel.length="+ sym_kernel.length);
System.out.println("asym_kernel.length="+asym_kernel.length); System.out.println("asym_kernel.length="+asym_kernel.length);
System.out.println("convolved.length="+convolved.length); System.out.println("convolved.length="+convolved.length);
}
SDFA_INSTANCE.showArrays(sym_kernel, DCT_PARAMETERS.dct_size, DCT_PARAMETERS.dct_size, "sym_kernel"); SDFA_INSTANCE.showArrays(sym_kernel, DCT_PARAMETERS.dct_size, DCT_PARAMETERS.dct_size, "sym_kernel");
SDFA_INSTANCE.showArrays(asym_kernel, DCT_PARAMETERS.asym_size, DCT_PARAMETERS.asym_size, "asym_kernel"); SDFA_INSTANCE.showArrays(asym_kernel, DCT_PARAMETERS.asym_size, DCT_PARAMETERS.asym_size, "asym_kernel");
SDFA_INSTANCE.showArrays(compare_kernels, target_expanded_size, target_expanded_size, true, "compare_kernels"); SDFA_INSTANCE.showArrays(compare_kernels, target_expanded_size, target_expanded_size, true, "compare_kernels");
...@@ -2938,13 +2963,29 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2938,13 +2963,29 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
if (!DCT_PARAMETERS.showDialog()) return; if (!DCT_PARAMETERS.showDialog()) return;
FactorConvKernel factorConvKernel = new FactorConvKernel(DCT_PARAMETERS.dbg_mode == 1); FactorConvKernel factorConvKernel = new FactorConvKernel(DCT_PARAMETERS.dbg_mode == 1);
factorConvKernel.setDebugLevel(DEBUG_LEVEL); factorConvKernel.setDebugLevel(DEBUG_LEVEL);
factorConvKernel.setTargetWindowMode(DCT_PARAMETERS.dbg_window_mode); factorConvKernel.setTargetWindowMode(DCT_PARAMETERS.dbg_window_mode, DCT_PARAMETERS.centerWindowToTarget);
factorConvKernel.numIterations = DCT_PARAMETERS.LMA_steps; factorConvKernel.numIterations = DCT_PARAMETERS.LMA_steps;
factorConvKernel.setAsymCompactness( factorConvKernel.setAsymCompactness(
DCT_PARAMETERS.compactness, DCT_PARAMETERS.compactness,
DCT_PARAMETERS.asym_tax_free); DCT_PARAMETERS.asym_tax_free);
int target_kernel_size = 2*DCT_PARAMETERS.dct_size - 1; int target_kernel_size = 2*DCT_PARAMETERS.dct_size - 1;
int target_expanded_size = 2*DCT_PARAMETERS.dct_size + DCT_PARAMETERS.asym_size -2;
double [] target_expanded = null;
if ((EYESIS_DCT != null) && EYESIS_DCT.kernelImageSet() && (DCT_PARAMETERS.color_channel >= 0)){
System.out.println("Using extracted target kernel");
double [] src_kernel = EYESIS_DCT.extractOneKernelFromStack(
CONVOLVE_FFT_SIZE/2 , // 64
DCT_PARAMETERS.color_channel, // 0..2
DCT_PARAMETERS.tileX, // horizontal number of kernel to extract
DCT_PARAMETERS.tileY); // vertical number of kernel to extract
target_expanded = EYESIS_DCT.reformatKernel(
src_kernel,// will be blured in-place
CONVOLVE_FFT_SIZE/2, // typical 64
target_expanded_size, // typical 15
DCT_PARAMETERS.decimation,// typical 2
DCT_PARAMETERS.decimateSigma);
} else {
System.out.println("Using synthesized target kernel");
double [] target_kernel = new double [target_kernel_size * target_kernel_size]; double [] target_kernel = new double [target_kernel_size * target_kernel_size];
for (int ii=0; ii < target_kernel.length; ii++) target_kernel[ii]=0.0; for (int ii=0; ii < target_kernel.length; ii++) target_kernel[ii]=0.0;
double dist = Math.sqrt((DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)*(DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)+ double dist = Math.sqrt((DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)*(DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)+
...@@ -2959,17 +3000,13 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2959,17 +3000,13 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
System.out.println(ii+": "+((DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)*ii/dist + DCT_PARAMETERS.dbg_x)+ System.out.println(ii+": "+((DCT_PARAMETERS.dbg_x1-DCT_PARAMETERS.dbg_x)*ii/dist + DCT_PARAMETERS.dbg_x)+
" / "+ ((DCT_PARAMETERS.dbg_y1-DCT_PARAMETERS.dbg_y)*ii/dist + DCT_PARAMETERS.dbg_y)+" ("+dbg_x+":"+dbg_y+")"); " / "+ ((DCT_PARAMETERS.dbg_y1-DCT_PARAMETERS.dbg_y)*ii/dist + DCT_PARAMETERS.dbg_y)+" ("+dbg_x+":"+dbg_y+")");
} }
} }
double blurSigma = DCT_PARAMETERS.dbg_sigma; double blurSigma = DCT_PARAMETERS.dbg_sigma;
DoubleGaussianBlur gb=null; DoubleGaussianBlur gb=null;
if (blurSigma>0) gb=new DoubleGaussianBlur(); if (blurSigma>0) gb=new DoubleGaussianBlur();
if (blurSigma>0) gb.blurDouble(target_kernel, target_kernel_size, target_kernel_size, blurSigma, blurSigma, 0.01); if (blurSigma>0) gb.blurDouble(target_kernel, target_kernel_size, target_kernel_size, blurSigma, blurSigma, 0.01);
// SDFA_INSTANCE.showArrays(target_kernel, target_kernel_size, target_kernel_size, "target_kernel"); target_expanded = new double [target_expanded_size * target_expanded_size];
int target_expanded_size = 2*DCT_PARAMETERS.dct_size + DCT_PARAMETERS.asym_size -2;
double [] target_expanded = new double [target_expanded_size * target_expanded_size];
for (int ii=0; ii < target_expanded.length; ii++) target_expanded[ii]=0.0; for (int ii=0; ii < target_expanded.length; ii++) target_expanded[ii]=0.0;
int left_top_margin = ((DCT_PARAMETERS.asym_size-1)/2); int left_top_margin = ((DCT_PARAMETERS.asym_size-1)/2);
for (int ii=0;ii < target_kernel_size; ii++){ for (int ii=0;ii < target_kernel_size; ii++){
...@@ -2978,7 +3015,8 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2978,7 +3015,8 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
target_kernel[ii*target_kernel_size + jj]; target_kernel[ii*target_kernel_size + jj];
} }
} }
}
long startTime = System.nanoTime();
int numPixels = factorConvKernel.calcKernels( int numPixels = factorConvKernel.calcKernels(
target_expanded, target_expanded,
DCT_PARAMETERS.asym_size, DCT_PARAMETERS.asym_size,
...@@ -2987,21 +3025,14 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -2987,21 +3025,14 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
DCT_PARAMETERS.asym_pixels, DCT_PARAMETERS.asym_pixels,
DCT_PARAMETERS.asym_distance, DCT_PARAMETERS.asym_distance,
DCT_PARAMETERS.seed_size); DCT_PARAMETERS.seed_size);
/* System.out.println(
public int calcKernels( "calcKernels() number of asym pixels = "+numPixels+
double []target_kernel, " RMS = "+factorConvKernel.getRMSes()[0]+
int asym_size, ", RMSPure = "+factorConvKernel.getRMSes()[1]+
int sym_radius, ", relRMSPure = "+(factorConvKernel.getRMSes()[1]/factorConvKernel.getTargetRMS())+
double fact_precision, ", number of LMA runs = "+factorConvKernel.getLMARuns()+
int asym_pixels, // maximal number of non-zero pixels in asymmmetrical kernel ", spent "+ IJ.d2s(0.000000001*(System.nanoTime()-startTime),3)+" sec");
int asym_distance){ // how far to seed a new pixel
*/
System.out.println("factorConvKernel.calcKernels() returned: >>> "+numPixels+ " <<<");
double [] sym_kernel = factorConvKernel.getSymKernel(); double [] sym_kernel = factorConvKernel.getSymKernel();
double [] asym_kernel = factorConvKernel.getAsymKernel(); double [] asym_kernel = factorConvKernel.getAsymKernel();
double [] convolved = factorConvKernel.getConvolved(); double [] convolved = factorConvKernel.getConvolved();
...@@ -3013,11 +3044,12 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -3013,11 +3044,12 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
weighted_diff100[ii] = diff100[ii]* target_weights[ii]; weighted_diff100[ii] = diff100[ii]* target_weights[ii];
} }
double [][] compare_kernels = {target_expanded, convolved, weighted_diff100,target_weights, diff100}; double [][] compare_kernels = {target_expanded, convolved, weighted_diff100,target_weights, diff100};
if (DEBUG_LEVEL>0) {
System.out.println("DCT_PARAMETERS.dct_size="+DCT_PARAMETERS.dct_size+" DCT_PARAMETERS.asym_size="+DCT_PARAMETERS.asym_size); System.out.println("DCT_PARAMETERS.dct_size="+DCT_PARAMETERS.dct_size+" DCT_PARAMETERS.asym_size="+DCT_PARAMETERS.asym_size);
System.out.println("sym_kernel.length="+ sym_kernel.length); System.out.println("sym_kernel.length="+ sym_kernel.length);
System.out.println("asym_kernel.length="+asym_kernel.length); System.out.println("asym_kernel.length="+asym_kernel.length);
System.out.println("convolved.length="+convolved.length); System.out.println("convolved.length="+convolved.length);
}
SDFA_INSTANCE.showArrays(sym_kernel, DCT_PARAMETERS.dct_size, DCT_PARAMETERS.dct_size, "sym_kernel"); SDFA_INSTANCE.showArrays(sym_kernel, DCT_PARAMETERS.dct_size, DCT_PARAMETERS.dct_size, "sym_kernel");
SDFA_INSTANCE.showArrays(asym_kernel, DCT_PARAMETERS.asym_size, DCT_PARAMETERS.asym_size, "asym_kernel"); SDFA_INSTANCE.showArrays(asym_kernel, DCT_PARAMETERS.asym_size, DCT_PARAMETERS.asym_size, "asym_kernel");
SDFA_INSTANCE.showArrays(compare_kernels, target_expanded_size, target_expanded_size, true, "compare_kernels"); SDFA_INSTANCE.showArrays(compare_kernels, target_expanded_size, target_expanded_size, true, "compare_kernels");
...@@ -3068,6 +3100,22 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -3068,6 +3100,22 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
EYESIS_CORRECTIONS.initSensorFiles(DEBUG_LEVEL); EYESIS_CORRECTIONS.initSensorFiles(DEBUG_LEVEL);
int numChannels=EYESIS_CORRECTIONS.getNumChannels(); int numChannels=EYESIS_CORRECTIONS.getNumChannels();
} else if (label.equals("Select kernels image")) {
if (!DCT_PARAMETERS.showDialog()) return;
if (EYESIS_DCT == null){
EYESIS_DCT = new EyesisDCT (
EYESIS_CORRECTIONS,
CORRECTION_PARAMETERS,
DCT_PARAMETERS);
}
ImagePlus imp_src = WindowManager.getCurrentImage();
if (imp_src==null){
IJ.showMessage("Error","3-layer files of deconvolution kernels is required");
return;
}
EYESIS_DCT.setKernelImageFile( imp_src);
} else if (label.equals("Create DCT kernels")) { } else if (label.equals("Create DCT kernels")) {
if (!DCT_PARAMETERS.showDialog()) return; if (!DCT_PARAMETERS.showDialog()) return;
if (EYESIS_DCT == null){ if (EYESIS_DCT == null){
...@@ -3105,18 +3153,14 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos ...@@ -3105,18 +3153,14 @@ private Panel panel1,panel2,panel3,panel4,panel5,panel5a, panel6,panel7,panelPos
EYESIS_CORRECTIONS.initSensorFiles(DEBUG_LEVEL); EYESIS_CORRECTIONS.initSensorFiles(DEBUG_LEVEL);
EYESIS_DCT.createDCTKernels(); EYESIS_DCT.createDCTKernels(
DCT_PARAMETERS,
CONVOLVE_FFT_SIZE/2,
/* THREADS_MAX,
EYESIS_CORRECTIONS.updateImageNoiseGains( UPDATE_STATUS, // update status info
NONLIN_PARAMETERS, //EyesisCorrectionParameters.NonlinParameters nonlinParameters, DEBUG_LEVEL);
CONVOLVE_FFT_SIZE, //int fftSize, // 128 - fft size, kernel size should be size/2 // EyesisCorrectionParameters.DCTParameters dCTParameters,
THREADS_MAX, // int threadsMax, // maximal number of threads to launch // int srcKernelSize,
UPDATE_STATUS, // boolean updateStatus,
DEBUG_LEVEL); //int globalDebugLevel){
*/
} }
DEBUG_LEVEL=MASTER_DEBUG_LEVEL; DEBUG_LEVEL=MASTER_DEBUG_LEVEL;
......
src/main/java/FactorConvKernel.java
View file @ 91616de2
...@@ -78,7 +78,7 @@ public class FactorConvKernel { ...@@ -78,7 +78,7 @@ public class FactorConvKernel {
public double currentRMS ; public double currentRMS ;
public double firstRMS; public double firstRMS;
public double nextRMS ; public double nextRMS ;
public double [] RMSes = {-1.0,-1.0};
public double currentRMSPure=-1.0; // calculated RMS for the currentVector->currentfX public double currentRMSPure=-1.0; // calculated RMS for the currentVector->currentfX
public double nextRMSPure= -1.0; // calculated RMS for the nextVector->nextfX public double nextRMSPure= -1.0; // calculated RMS for the nextVector->nextfX
public double firstRMSPure= -1.0; // RMS before current series of LMA started public double firstRMSPure= -1.0; // RMS before current series of LMA started
...@@ -98,7 +98,8 @@ public class FactorConvKernel { ...@@ -98,7 +98,8 @@ public class FactorConvKernel {
public double goal_rms_pure; public double goal_rms_pure;
public LMAData lMAData = null; public LMAData lMAData = null;
public int numLMARuns = 0; public int numLMARuns = 0;
public int target_window_mode = 2; // 0 - none, 1 - square, 2 - sin public int target_window_mode = 2; // 0 - none, 1 - square, 2 - sin, 3 - sin^2
public boolean centerWindowToTarget = true; // center convolution weights around target kernel center
public class LMAArrays { public class LMAArrays {
public double [][] jTByJ= null; // jacobian multiplied by Jacobian transposed public double [][] jTByJ= null; // jacobian multiplied by Jacobian transposed
...@@ -141,6 +142,7 @@ public class FactorConvKernel { ...@@ -141,6 +142,7 @@ public class FactorConvKernel {
public double [] saved_rmses = null; public double [] saved_rmses = null;
public double [] first_rmses = null; public double [] first_rmses = null;
public int target_window_mode = 2; // 0 - none, 1 - square, 2 - sin public int target_window_mode = 2; // 0 - none, 1 - square, 2 - sin
public boolean centerWindowToTarget = true; // center convolution weights around target kernel center
public LMAData(){ public LMAData(){
...@@ -148,8 +150,9 @@ public class FactorConvKernel { ...@@ -148,8 +150,9 @@ public class FactorConvKernel {
public LMAData( int debugLevel){ public LMAData( int debugLevel){
this.debugLevel = debugLevel; this.debugLevel = debugLevel;
} }
public void setTargetWindowMode(int mode){ public void setTargetWindowMode(int mode, boolean centerWindowToTarget){
target_window_mode = mode; this.target_window_mode = mode;
this.centerWindowToTarget = centerWindowToTarget;
} }
public void initLMAData(){ public void initLMAData(){
...@@ -390,20 +393,36 @@ public class FactorConvKernel { ...@@ -390,20 +393,36 @@ public class FactorConvKernel {
System.out.println("getWeight(), delete jacobian"); System.out.println("getWeight(), delete jacobian");
} }
if ((weights[0] == null) || (weights[0].length != target_kernel.length)){ if ((weights[0] == null) || (weights[0].length != target_kernel.length)){
if (debugLevel > 2){ int xc = 0;
System.out.println("Convolution weight is null/outdated, regenerating default (only center part)"); int yc = 0;
int conv_size = asym_size + 2*sym_radius-2;
int cc = conv_size/2;
if (this.centerWindowToTarget) {
double s0=0.0,sx=0.0,sy=0.0;
for (int i = 0; i < conv_size; i++){
for (int j = 0; j < conv_size; j++){
double d = target_kernel[conv_size*i+j];
d *= d;
s0 += d;
sx += d* (j - cc);
sy += d* (i - cc);
}
}
xc = (int) Math.round(sx/s0);
yc = (int) Math.round(sy/s0);
} }
double [] sins = new double [2*sym_radius-1]; double [] sins = new double [2*sym_radius-1];
for (int i = 1; i< 2*sym_radius; i++) sins[i-1] = Math.sin(Math.PI*i/(2.0 * sym_radius)); for (int i = 1; i< 2*sym_radius; i++) sins[i-1] = Math.sin(Math.PI*i/(2.0 * sym_radius));
weights[0] = new double [target_kernel.length]; weights[0] = new double [target_kernel.length];
int conv_size = asym_size + 2*sym_radius-2; int left_margin = ((asym_size-1)/2) +xc; // inclusive
int left_top_margin = ((asym_size-1)/2); // inclusive int top_margin = ((asym_size-1)/2) + yc; // inclusive
int right_bottom_margin = left_top_margin + (2 * sym_radius - 1); // exclusive int right_margin = left_margin + (2 * sym_radius - 1); // exclusive
int bottom_margin = top_margin + (2 * sym_radius - 1); // exclusive
for (int i = 0; i < conv_size; i++) { for (int i = 0; i < conv_size; i++) {
for (int j = 0; j < conv_size; j++){ for (int j = 0; j < conv_size; j++){
int cindx = i * conv_size + j; int cindx = i * conv_size + j;
if ((i >= left_top_margin) && (i < right_bottom_margin) && (j >= left_top_margin) && (j < right_bottom_margin)) { if ((i >= top_margin) && (i < bottom_margin) && (j >= left_margin) && (j < right_margin)) {
weights[0][cindx] = (target_window_mode>=2)?(sins[i-left_top_margin]*sins[j-left_top_margin]):1.0; weights[0][cindx] = (target_window_mode>=2)?(sins[i-top_margin]*sins[j-left_margin]):1.0;
if (target_window_mode == 3) { if (target_window_mode == 3) {
weights[0][cindx]*=weights[0][cindx]; weights[0][cindx]*=weights[0][cindx];
} }
...@@ -413,6 +432,8 @@ public class FactorConvKernel { ...@@ -413,6 +432,8 @@ public class FactorConvKernel {
} }
} }
} }
// public boolean centerWindowToTarget = true; // center convolution weights around target kernel center
//target_window_mode //target_window_mode
rebuildMapsPars(false); // only if it does not exist rebuildMapsPars(false); // only if it does not exist
double sw = 0.0; double sw = 0.0;
...@@ -749,6 +770,7 @@ public class FactorConvKernel { ...@@ -749,6 +770,7 @@ public class FactorConvKernel {
} // class LMAData } // class LMAData
public FactorConvKernel(){ public FactorConvKernel(){
this.new_mode =true;
} }
public FactorConvKernel(boolean new_mode){ public FactorConvKernel(boolean new_mode){
...@@ -760,8 +782,17 @@ public class FactorConvKernel { ...@@ -760,8 +782,17 @@ public class FactorConvKernel {
this.sym_radius = sym_radius; this.sym_radius = sym_radius;
} }
public void setTargetWindowMode(int mode){ public void setTargetWindowMode(int mode, boolean centerWindowToTarget){
target_window_mode = mode; target_window_mode = mode;
this.centerWindowToTarget = centerWindowToTarget;
}
public double[] getRMSes(){
return this.RMSes; // lMAData.getRMSes();
}
public double getTargetRMS(){
return lMAData.getTargetRMSW();
} }
public double [] generateAsymWeights( public double [] generateAsymWeights(
...@@ -844,7 +875,9 @@ public class FactorConvKernel { ...@@ -844,7 +875,9 @@ public class FactorConvKernel {
lMAData.setSymKernel(null, sym_mask_frozen); lMAData.setSymKernel(null, sym_mask_frozen);
levenbergMarquardt(); levenbergMarquardt();
lMAData.setSymKernel(null, sym_mask); lMAData.setSymKernel(null, sym_mask);
return levenbergMarquardt(); boolean OK = levenbergMarquardt();
this.RMSes = lMAData.getRMSes().clone();
return OK;
} else{ } else{
initLevenbergMarquardt_old(fact_precision, seed_size, asym_random); initLevenbergMarquardt_old(fact_precision, seed_size, asym_random);
if (mask != null){ if (mask != null){
...@@ -1116,6 +1149,7 @@ public class FactorConvKernel { ...@@ -1116,6 +1149,7 @@ public class FactorConvKernel {
System.out.println(); System.out.println();
} }
} }
this.RMSes = RMSes.clone();
return numAsym; return numAsym;
} }
...@@ -2148,16 +2182,7 @@ public class FactorConvKernel { ...@@ -2148,16 +2182,7 @@ public class FactorConvKernel {
private void initLevenbergMarquardt(double fact_precision, int seed_size, double asym_random){ private void initLevenbergMarquardt(double fact_precision, int seed_size, double asym_random){
lMAData = new LMAData(debugLevel); lMAData = new LMAData(debugLevel);
lMAData.setTarget(target_kernel); lMAData.setTarget(target_kernel);
lMAData.setTargetWindowMode(target_window_mode); lMAData.setTargetWindowMode(target_window_mode, centerWindowToTarget);
/* double s = 0.0;
for (int i = 0; i<target_kernel.length; i++){
s+= target_kernel[i]*target_kernel[i];
}
this.goal_rms_pure = Math.sqrt(s/target_kernel.length)*fact_precision;
*/
// this.currentVector = setInitialVector(target_kernel, null); // should be (asym_size + 2*sym_radius-1)**2
double [][]kernels = setInitialVector(target_kernel, seed_size, asym_random); // should be (asym_size + 2*sym_radius-1)**2 double [][]kernels = setInitialVector(target_kernel, seed_size, asym_random); // should be (asym_size + 2*sym_radius-1)**2
sym_kernel_scale = kernels[0][0]; sym_kernel_scale = kernels[0][0];
for (int i=0;i<kernels[0].length;i++) if (!Double.isNaN(kernels[0][i])) kernels[0][i] /=sym_kernel_scale; for (int i=0;i<kernels[0].length;i++) if (!Double.isNaN(kernels[0][i])) kernels[0][i] /=sym_kernel_scale;
...@@ -2166,13 +2191,9 @@ public class FactorConvKernel { ...@@ -2166,13 +2191,9 @@ public class FactorConvKernel {
lMAData.setAsymKernel(kernels[1]); lMAData.setAsymKernel(kernels[1]);
this.goal_rms_pure = lMAData.getTargetRMSW()*fact_precision; this.goal_rms_pure = lMAData.getTargetRMSW()*fact_precision;
this.target_rms = lMAData.getTargetRMSW(); //Math.sqrt(s/target_kernel.length); this.target_rms = lMAData.getTargetRMSW(); //Math.sqrt(s/target_kernel.length);
// lMAData.invalidateLMAArrays(); // should be in each run , not at init
resetLMARuns(); resetLMARuns();
// this.currentVector = setVectorFromKernels(kernels[0], kernels[1]);
} }
private boolean levenbergMarquardt(){ private boolean levenbergMarquardt(){
long startTime=System.nanoTime(); long startTime=System.nanoTime();
this.iterationStepNumber = 0; this.iterationStepNumber = 0;
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment