Skip to content
Projects
Groups
Snippets
Help
Loading...
Help
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
I
imagej-elphel
Project
Project
Details
Activity
Releases
Cycle Analytics
Repository
Repository
Files
Commits
Branches
Tags
Contributors
Graph
Compare
Charts
Issues
3
Issues
3
List
Board
Labels
Milestones
Wiki
Wiki
Members
Members
Collapse sidebar
Close sidebar
Activity
Graph
Charts
Create a new issue
Commits
Issue Boards
Open sidebar
Elphel
imagej-elphel
Commits
3479dfe6
Commit
3479dfe6
authored
Jul 09, 2014
by
Andrey Filippov
Browse files
Options
Browse Files
Download
Email Patches
Plain Diff
working on saving/restoring new parameters
parent
ca5bb1b6
Changes
1
Show whitespace changes
Inline
Side-by-side
Showing
1 changed file
with
3813 additions
and
3650 deletions
+3813
-3650
FocusingField.java
src/main/java/FocusingField.java
+3813
-3650
No files found.
src/main/java/FocusingField.java
View file @
3479dfe6
/**
**
** FocusingField.jave - save/restore/process sagittal/tangential PSF width
** over FOV, together with related data
**
** Copyright (C) 2014 Elphel, Inc.
**
** -----------------------------------------------------------------------------**
**
**
CalibrationHardwareInterface.java is free software: you can redistribute it and/or modify
**
it under the terms of the GNU General Public License as published by
**
the Free Software Foundation, either version 3 of the License, or
**
(at your option) any later version.
**
**
This program is distributed in the hope that it will be useful,
**
but WITHOUT ANY WARRANTY; without even the implied warranty of
** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
See the
**
GNU General Public License for more details.
**
**
You should have received a copy of the GNU General Public License
** along with this program.
If not, see <http://www.gnu.org/licenses/>.
** -----------------------------------------------------------------------------**
**
*/
**
** FocusingField.jave - save/restore/process sagittal/tangential PSF width
** over FOV, together with related data
**
** Copyright (C) 2014 Elphel, Inc.
**
** -----------------------------------------------------------------------------**
**
**
CalibrationHardwareInterface.java is free software: you can redistribute it and/or modify
**
it under the terms of the GNU General Public License as published by
**
the Free Software Foundation, either version 3 of the License, or
**
(at your option) any later version.
**
**
This program is distributed in the hope that it will be useful,
**
but WITHOUT ANY WARRANTY; without even the implied warranty of
** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
See the
**
GNU General Public License for more details.
**
**
You should have received a copy of the GNU General Public License
** along with this program.
If not, see <http://www.gnu.org/licenses/>.
** -----------------------------------------------------------------------------**
**
*/
import
ij.IJ
;
import
ij.gui.GenericDialog
;
...
...
@@ -33,6 +33,7 @@ import java.io.FileWriter;
import
java.io.IOException
;
import
java.util.ArrayList
;
import
java.util.Arrays
;
import
java.util.Properties
;
import
java.util.concurrent.atomic.AtomicInteger
;
import
org.apache.commons.configuration.ConfigurationException
;
...
...
@@ -44,6 +45,7 @@ import org.apache.commons.configuration.XMLConfiguration;
//import Distortions.LMAArrays; // may still reuse?
import
Jama.LUDecomposition
;
import
Jama.Matrix
;
...
...
@@ -83,91 +85,92 @@ public class FocusingField {
double
[]
currentVector
=
null
;
double
[]
nextVector
=
null
;
double
[]
savedVector
=
null
;
private
LMAArrays
lMAArrays
=
null
;
private
LMAArrays
savedLMAArrays
=
null
;
private
double
[]
currentfX
=
null
;
// array of "f(x)" - simulated data for all images, combining pixel-X and pixel-Y (odd/even)
private
double
[]
nextfX
=
null
;
// array of "f(x)" - simulated data for all images, combining pixel-X and pixel-Y (odd/even)
private
double
currentRMS
=-
1.0
;
// calculated RMS for the currentVector->currentfX
private
double
currentRMSPure
=-
1.0
;
// calculated RMS for the currentVector->currentfX
private
double
nextRMS
=-
1.0
;
// calculated RMS for the nextVector->nextfX
private
double
nextRMSPure
=-
1.0
;
// calculated RMS for the nextVector->nextfX
private
double
firstRMS
=-
1.0
;
// RMS before current series of LMA started
private
double
firstRMSPure
=-
1.0
;
// RMS before current series of LMA started
private
double
lambdaStepUp
=
8.0
;
// multiply lambda by this if result is worse
private
double
lambdaStepDown
=
0.5
;
// multiply lambda by this if result is better
private
double
thresholdFinish
=
0.001
;
// (copied from series) stop iterations if 2 last steps had less improvement (but not worsening )
private
int
numIterations
=
100
;
// maximal number of iterations
private
double
maxLambda
=
100.0
;
// max lambda to fail
private
double
lambda
=
0.001
;
// copied from series
private
double
[]
lastImprovements
=
{-
1.0
,-
1.0
};
// {last improvement, previous improvement}. If both >0 and < thresholdFinish - done
private
int
iterationStepNumber
=
0
;
private
boolean
stopEachStep
=
true
;
// open dialog after each fitting step
private
boolean
stopOnFailure
=
true
;
// open dialog when fitting series failed
private
boolean
showParams
=
false
;
// show modified parameters
private
boolean
showDisabledParams
=
false
;
private
boolean
showCorrectionParams
=
false
;
private
long
startTime
=
0
;
private
AtomicInteger
stopRequested
=
null
;
// 1 - stop now, 2 - when convenient
private
boolean
saveSeries
=
false
;
// just for the dialog
private
boolean
showMotors
=
true
;
private
boolean
[]
showMeasCalc
=
{
true
,
true
,
true
};
private
boolean
[]
showColors
=
{
true
,
true
,
true
};
private
boolean
[]
showDirs
=
{
true
,
true
};
private
boolean
[]
showSamples
=
null
;
private
boolean
showAllSamples
=
true
;
private
boolean
showIgnoredData
=
false
;
private
boolean
showRad
=
true
;
private
boolean
[][][][][]
sampleMask
=
null
;
private
boolean
correct_measurement_ST
=
true
;
private
boolean
updateWeightWhileFitting
=
false
;
public
int
debugLevel
;
public
boolean
debugDerivatives
;
public
boolean
debugDerivativesFxDxDy
=
false
;
private
double
k_red
=
0.7
;
private
double
k_blue
=
0.4
;
private
double
qb_scan_below
=-
20.0
;
// um
private
double
qb_scan_above
=
60.0
;
// um
private
double
qb_scan_step
=
0.5
;
// um
private
boolean
qb_use_corrected
=
true
;
private
boolean
qb_invert
=
true
;
private
boolean
rslt_show_z_axial
=
true
;
private
boolean
rslt_show_z_individual
=
true
;
private
boolean
rslt_show_f_axial
=
true
;
private
boolean
rslt_show_f_individual
=
true
;
private
double
rslt_scan_below
=-
10.0
;
private
double
rslt_scan_above
=
10.0
;
private
double
rslt_scan_step
=
5.0
;
private
boolean
rslt_mtf50_mode
=
true
;
private
boolean
[]
rslt_show_chn
={
true
,
true
,
true
,
true
,
true
,
true
};
private
LMAArrays
lMAArrays
=
null
;
private
LMAArrays
savedLMAArrays
=
null
;
private
double
[]
currentfX
=
null
;
// array of "f(x)" - simulated data for all images, combining pixel-X and pixel-Y (odd/even)
private
double
[]
nextfX
=
null
;
// array of "f(x)" - simulated data for all images, combining pixel-X and pixel-Y (odd/even)
private
double
currentRMS
=-
1.0
;
// calculated RMS for the currentVector->currentfX
private
double
currentRMSPure
=-
1.0
;
// calculated RMS for the currentVector->currentfX
private
double
nextRMS
=-
1.0
;
// calculated RMS for the nextVector->nextfX
private
double
nextRMSPure
=-
1.0
;
// calculated RMS for the nextVector->nextfX
private
double
firstRMS
=-
1.0
;
// RMS before current series of LMA started
private
double
firstRMSPure
=-
1.0
;
// RMS before current series of LMA started
private
double
lambdaStepUp
=
8.0
;
// multiply lambda by this if result is worse
private
double
lambdaStepDown
=
0.5
;
// multiply lambda by this if result is better
private
double
thresholdFinish
=
0.001
;
// (copied from series) stop iterations if 2 last steps had less improvement (but not worsening )
private
int
numIterations
=
100
;
// maximal number of iterations
private
double
maxLambda
=
100.0
;
// max lambda to fail
private
double
lambda
=
0.001
;
// copied from series
private
double
[]
lastImprovements
=
{-
1.0
,-
1.0
};
// {last improvement, previous improvement}. If both >0 and < thresholdFinish - done
private
int
iterationStepNumber
=
0
;
private
boolean
stopEachStep
=
true
;
// open dialog after each fitting step
private
boolean
stopOnFailure
=
true
;
// open dialog when fitting series failed
private
boolean
showParams
=
false
;
// show modified parameters
private
boolean
showDisabledParams
=
false
;
private
boolean
showCorrectionParams
=
false
;
private
boolean
keepCorrectionParameters
=
true
;
private
long
startTime
=
0
;
private
AtomicInteger
stopRequested
=
null
;
// 1 - stop now, 2 - when convenient
private
boolean
saveSeries
=
false
;
// just for the dialog
private
boolean
showMotors
=
true
;
private
boolean
[]
showMeasCalc
=
{
true
,
true
,
true
};
private
boolean
[]
showColors
=
{
true
,
true
,
true
};
private
boolean
[]
showDirs
=
{
true
,
true
};
private
boolean
[]
showSamples
=
null
;
private
boolean
showAllSamples
=
true
;
private
boolean
showIgnoredData
=
false
;
private
boolean
showRad
=
true
;
private
boolean
[][][][][]
sampleMask
=
null
;
private
boolean
correct_measurement_ST
=
true
;
private
boolean
updateWeightWhileFitting
=
false
;
public
int
debugLevel
;
public
boolean
debugDerivatives
;
public
boolean
debugDerivativesFxDxDy
=
false
;
private
double
k_red
=
0.7
;
private
double
k_blue
=
0.4
;
private
double
qb_scan_below
=-
20.0
;
// um
private
double
qb_scan_above
=
60.0
;
// um
private
double
qb_scan_step
=
0.5
;
// um
private
boolean
qb_use_corrected
=
true
;
private
boolean
qb_invert
=
true
;
private
boolean
rslt_show_z_axial
=
true
;
private
boolean
rslt_show_z_individual
=
true
;
private
boolean
rslt_show_f_axial
=
true
;
private
boolean
rslt_show_f_individual
=
true
;
private
double
rslt_scan_below
=-
10.0
;
private
double
rslt_scan_above
=
10.0
;
private
double
rslt_scan_step
=
5.0
;
private
boolean
rslt_mtf50_mode
=
true
;
private
boolean
[]
rslt_show_chn
={
true
,
true
,
true
,
true
,
true
,
true
};
// public double fwhm_to_mtf50=500.0; // put actual number
public
double
fwhm_to_mtf50
=
2
*
Math
.
log
(
2.0
)/
Math
.
PI
*
1000
;
//pi/0.004
public
double
fwhm_to_mtf50
=
2
*
Math
.
log
(
2.0
)/
Math
.
PI
*
1000
;
//pi/0.004
public
void
setDebugLevel
(
int
debugLevel
){
this
.
debugLevel
=
debugLevel
;
}
public
class
LMAArrays
{
// reuse from Distortions?
public
double
[][]
jTByJ
=
null
;
// jacobian multiplied by Jacobian transposed
public
double
[]
jTByDiff
=
null
;
// jacobian multiplied difference vector
public
LMAArrays
clone
()
{
public
class
LMAArrays
{
// reuse from Distortions?
public
double
[][]
jTByJ
=
null
;
// jacobian multiplied by Jacobian transposed
public
double
[]
jTByDiff
=
null
;
// jacobian multiplied difference vector
public
LMAArrays
clone
()
{
LMAArrays
lma
=
new
LMAArrays
();
lma
.
jTByJ
=
this
.
jTByJ
.
clone
();
for
(
int
i
=
0
;
i
<
this
.
jTByJ
.
length
;
i
++)
lma
.
jTByJ
[
i
]=
this
.
jTByJ
[
i
].
clone
();
lma
.
jTByDiff
=
this
.
jTByDiff
.
clone
();
return
lma
;
}
}
}
}
public
enum
MECH_PAR
{
public
enum
MECH_PAR
{
K0
,
// Average motor center travel","um/step","0.0124"},
KD1
,
// M1 and M2 travel disbalance","um/step","0.0"},
KD3
,
// M3 to average of M1 and M2 travel disbalance","um/step","0.0"},
...
...
@@ -184,20 +187,20 @@ public class FocusingField {
z0
,
// center shift, positive away form the lens","um","0.0"},
tx
,
// horizontal tilt", "um/pix","0.0"},
ty
// vertical tilt", "um/pix","0.0"}};
}
public
static
double
getPixelMM
(){
return
PIXEL_SIZE
;}
public
static
double
getPixelUM
(){
return
PIXEL_SIZE
*
1000
;}
public
int
flattenIndex
(
int
i
,
int
j
){
return
j
+
i
*
sampleCoord
[
0
].
length
;}
public
int
getNumSamples
(){
return
sampleCoord
.
length
*
sampleCoord
[
0
].
length
;}
public
int
getSampleWidth
(){
return
sampleCoord
[
0
].
length
;}
public
double
[][]
flattenSampleCoord
(){
}
public
static
double
getPixelMM
(){
return
PIXEL_SIZE
;}
public
static
double
getPixelUM
(){
return
PIXEL_SIZE
*
1000
;}
public
int
flattenIndex
(
int
i
,
int
j
){
return
j
+
i
*
sampleCoord
[
0
].
length
;}
public
int
getNumSamples
(){
return
sampleCoord
.
length
*
sampleCoord
[
0
].
length
;}
public
int
getSampleWidth
(){
return
sampleCoord
[
0
].
length
;}
public
double
[][]
flattenSampleCoord
(){
///sampleCoord
double
[][]
flatSampleCoord
=
new
double
[
sampleCoord
.
length
*
sampleCoord
[
0
].
length
][
2
];
int
index
=
0
;
for
(
int
i
=
0
;
i
<
sampleCoord
.
length
;
i
++)
for
(
int
j
=
0
;
j
<
sampleCoord
[
0
].
length
;
j
++)
flatSampleCoord
[
index
++]=
sampleCoord
[
i
][
j
];
return
flatSampleCoord
;
// last dimension is not cloned
}
public
class
MeasuredSample
{
}
public
class
MeasuredSample
{
public
int
[]
motors
=
new
int
[
3
];
public
String
timestamp
;
public
double
px
;
...
...
@@ -230,8 +233,8 @@ public class FocusingField {
this
.
dPxyc
[
1
]=
dPyc
;
this
.
sampleIndex
=
sampleIndex
;
}
}
public
boolean
configureDataVector
(
String
title
,
boolean
forcenew
){
}
public
boolean
configureDataVector
(
String
title
,
boolean
forcenew
){
if
((
fieldFitting
==
null
)
&&
!
forcenew
){
forcenew
=
true
;
}
...
...
@@ -313,13 +316,13 @@ public class FocusingField {
this
.
savedVector
=
fieldFitting
.
createParameterVector
(
sagittalMaster
);
// initialVector
return
true
;
}
}
// includes deselected channels
public
void
setDataVector
(
MeasuredSample
[]
vector
){
// remove unused channels if any
// includes deselected channels
public
void
setDataVector
(
MeasuredSample
[]
vector
){
// remove unused channels if any
if
(
debugLevel
>
1
)
System
.
out
.
println
(
"+++++ (Re)calculating sample weights +++++"
);
boolean
[]
chanSel
=
fieldFitting
.
getSelectedChannels
();
int
numSamples
=
0
;
...
...
@@ -330,7 +333,7 @@ public class FocusingField {
int
corrLength
=
fieldFitting
.
getNumberOfCorrParameters
();
dataValues
=
new
double
[
dataVector
.
length
+
corrLength
];
dataWeights
=
new
double
[
dataVector
.
length
+
corrLength
];
//
sumWeights=0.0;
//
sumWeights=0.0;
int
mode
=
weighMode
;
double
kw
=
(
weightRadius
>
0.0
)?(-
0.5
*
getPixelMM
()*
getPixelMM
()/(
weightRadius
*
weightRadius
)):
0
;
//weightRadius
...
...
@@ -356,11 +359,11 @@ public class FocusingField {
dataValues
[
i
+
dataVector
.
length
]=
0.0
;
// correction target is always 0
dataWeights
[
i
+
dataVector
.
length
]=
1.0
;
// improve?
}
}
}
// for compatibility with Distortions class\
public
void
commitParameterVector
(
double
[]
vector
){
public
void
commitParameterVector
(
double
[]
vector
){
fieldFitting
.
commitParameterVector
(
vector
,
sagittalMaster
);
// recalculate measured S,T (depend on center) if center is among fitted parameters
boolean
[]
centerSelect
=
fieldFitting
.
getCenterSelect
();
...
...
@@ -378,15 +381,15 @@ public class FocusingField {
}
}
}
}
public
double
[]
createFXandJacobian
(
double
[]
vector
,
boolean
createJacobian
){
public
double
[]
createFXandJacobian
(
double
[]
vector
,
boolean
createJacobian
){
commitParameterVector
(
vector
);
return
createFXandJacobian
(
createJacobian
);
}
}
public
double
[]
createFXandJacobian
(
boolean
createJacobian
){
public
double
[]
createFXandJacobian
(
boolean
createJacobian
){
int
numCorrPar
=
fieldFitting
.
getNumberOfCorrParameters
();
double
[]
fx
=
new
double
[
dataVector
.
length
+
numCorrPar
];
double
[][]
derivs
=
null
;
...
...
@@ -431,15 +434,6 @@ public class FocusingField {
// jacobian[i][n]=derivs[selChanIndices[ms.channel]][i];
jacobian
[
i
][
n
]=
thisDerivs
[
i
];
}
/*
if ((fieldFitting.sampleCorrChnParIndex!=null) && (fieldFitting.sampleCorrChnParIndex[ms.channel]!=null)){
for (int i=0;i<fieldFitting.getNumCurvars()[0];i++){
if (fieldFitting.sampleCorrChnParIndex[ms.channel][i]>=0){
jacobian[i][n]
}
}
}
*/
//TODO: correct /dpx, /dpy to compensate for measured S,T calcualtion
boolean
[]
centerSelect
=
fieldFitting
.
getCenterSelect
();
if
(
correct_measurement_ST
&&
(
centerSelect
[
0
]
||
centerSelect
[
1
])){
// do not do that if both X and Y are disabled
...
...
@@ -460,10 +454,10 @@ public class FocusingField {
for
(
int
i
=
0
;
i
<
numSamples
;
i
++){
double
f
=
0.0
;
for
(
int
j
=
0
;
j
<
numSamples
;
j
++){
f
+=
fieldFitting
.
sampleCor
r
[
pindex
+
j
]*
fieldFitting
.
sampleCorrCrossWeights
[
chn
][
np
][
i
][
j
];
f
+=
fieldFitting
.
sampleCorrVecto
r
[
pindex
+
j
]*
fieldFitting
.
sampleCorrCrossWeights
[
chn
][
np
][
i
][
j
];
}
fx
[
index
]=
f
;
// f+=fieldFitting.sampleCor
r[pindex+i]
// f+=fieldFitting.sampleCorrVecto
r[pindex+i]
if
(
createJacobian
)
{
for
(
int
j
=
0
;
j
<
numSamples
;
j
++){
jacobian
[
numRegPars
+
pindex
+
j
][
index
]=
fieldFitting
.
sampleCorrCrossWeights
[
chn
][
np
][
i
][
j
];
...
...
@@ -477,9 +471,9 @@ public class FocusingField {
}
}
return
fx
;
}
}
public
double
getRMS
(
double
[]
fx
,
boolean
pure
){
public
double
getRMS
(
double
[]
fx
,
boolean
pure
){
int
len
=
pure
?
dataVector
.
length
:
fx
.
length
;
double
sum
=
0.0
;
double
sum_w
=
0.0
;
...
...
@@ -500,40 +494,40 @@ public class FocusingField {
sum
/=
sum_w
;
}
return
Math
.
sqrt
(
sum
);
}
}
public
MeasuredSample
[]
createDataVector
(){
public
MeasuredSample
[]
createDataVector
(){
return
createDataVector
(
true
,
// boolean updateSelection,
currentPX0
,
// double centerPX,
currentPY0
);
//double centerPY
}
public
MeasuredSample
[]
createDataVector
(
}
public
MeasuredSample
[]
createDataVector
(
boolean
updateSelection
,
double
centerPX
,
double
centerPY
){
// use this data
return
createDataVector
(
return
createDataVector
(
updateSelection
,
centerPX
,
centerPY
,
(
this
.
useMinMeas
?
this
.
minMeas
:
null
),
// pixels
(
this
.
useMaxMeas
?
this
.
maxMeas
:
null
),
// pixels
(
this
.
useThresholdMax
?
this
.
thresholdMax
:
null
));
// pixels
}
}
/**
* Generate of usable measurement samples
* @param minMeas minimal measurement PSF radius (in pixels) - correction that increases result
* resolution in "sharp" areas (to compensate for measurement error). Individual for each of
*
the 6 color/direction components.
* @param updateSelection when true - updates selection of "good" samples, when false - reuses existent one
* @param maxMeas maximal measurement PSF radius (in pixels) - correction that decreases result
* resolution out-of-focus areas to compensate for the limited size of the PSF window.
* Individual for each of the 6 color/direction components.
* @param thresholdMax maximal PSF radius to consider data usable
* @return array of the MeasuredSample instances, including motors, PSF radius, channel and value
*/
/**
* Generate of usable measurement samples
* @param minMeas minimal measurement PSF radius (in pixels) - correction that increases result
* resolution in "sharp" areas (to compensate for measurement error). Individual for each of
*
the 6 color/direction components.
* @param updateSelection when true - updates selection of "good" samples, when false - reuses existent one
* @param maxMeas maximal measurement PSF radius (in pixels) - correction that decreases result
* resolution out-of-focus areas to compensate for the limited size of the PSF window.
* Individual for each of the 6 color/direction components.
* @param thresholdMax maximal PSF radius to consider data usable
* @return array of the MeasuredSample instances, including motors, PSF radius, channel and value
*/
/*
Need to find partial derivatives of each of the 3 coefficients: c2, s2 and cs for both px0 and py0
r2 = (x-x0)^2+ (y-y0)^2
...
...
@@ -569,8 +563,8 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
d_cs/dx0= delta_y*(2*delta_x^2-r2)/r2^2
d_cs/dy0= delta_x*(2*delta_y^2-r2)/r2^2
*/
public
MeasuredSample
[]
createDataVector
(
*/
public
MeasuredSample
[]
createDataVector
(
boolean
updateSelection
,
double
centerPX
,
double
centerPY
,
...
...
@@ -602,7 +596,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
2*d_cs/dx0= 2*delta_y*(2*delta_x^2-r2)/r2^2
2*d_cs/dy0= 2*delta_x*(2*delta_y^2-r2)/r2^2
*/
*/
double
[][][]
cosSin2Tab
=
new
double
[
sampleCoord
.
length
][
sampleCoord
[
0
].
length
][
3
];
double
[][][]
debugCosSin2Tab_dx
=
null
;
double
[][][]
debugCosSin2Tab_dy
=
null
;
...
...
@@ -883,7 +877,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
if
(
debugLevel
>
3
)
System
.
out
.
println
();
if
(
debugLevel
>
0
)
System
.
out
.
println
(
"createDataVector -> "
+
sampleList
.
size
()+
" elements"
);
return
sampleList
.
toArray
(
new
MeasuredSample
[
0
]);
}
}
/**
* Calculate differences vector
* @param fX vector of calculated pixelX,pixelY on the sensors
...
...
@@ -1086,12 +1080,12 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
// this.savedJacobian=this.jacobian;
this
.
savedLMAArrays
=
lMAArrays
.
clone
();
this
.
jacobian
=
null
;
// not needed, just to catch bugs
this
.
nextfX
=
createFXandJacobian
(
this
.
nextVector
,
true
);
this
.
lMAArrays
=
calculateJacobianArrays
(
this
.
nextfX
);
this
.
nextRMS
=
calcErrorDiffY
(
this
.
nextfX
,
false
);
this
.
nextRMSPure
=
calcErrorDiffY
(
this
.
nextfX
,
true
);
this
.
savedLMAArrays
=
lMAArrays
.
clone
();
this
.
jacobian
=
null
;
// not needed, just to catch bugs
this
.
nextfX
=
createFXandJacobian
(
this
.
nextVector
,
true
);
this
.
lMAArrays
=
calculateJacobianArrays
(
this
.
nextfX
);
this
.
nextRMS
=
calcErrorDiffY
(
this
.
nextfX
,
false
);
this
.
nextRMSPure
=
calcErrorDiffY
(
this
.
nextfX
,
true
);
this
.
lastImprovements
[
1
]=
this
.
lastImprovements
[
0
];
this
.
lastImprovements
[
0
]=
this
.
currentRMS
-
this
.
nextRMS
;
if
(
debugLevel
>
2
)
{
...
...
@@ -1142,9 +1136,9 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
System
.
out
.
println
(
"stepLevenbergMarquardtFirst("
+
debugLevel
+
")=>"
+
status
[
0
]+
","
+
status
[
1
]);
}
return
status
;
}
}
public
void
stepLevenbergMarquardtAction
(
int
debugLevel
){
//
public
void
stepLevenbergMarquardtAction
(
int
debugLevel
){
//
this
.
iterationStepNumber
++;
// apply/revert,modify lambda
if
(
debugLevel
>
1
)
{
...
...
@@ -1165,16 +1159,18 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
this
.
lambda
*=
this
.
lambdaStepUp
;
this
.
lMAArrays
=
this
.
savedLMAArrays
;
// restore Jt*J and Jt*diff
}
}
}
/**
* Dialog to select Levenberg-Marquardt algorithm and related parameters
* @return true if OK, false if canceled
*/
public
boolean
selectLMAParameters
(){
/**
* Dialog to select Levenberg-Marquardt algorithm and related parameters
* @return true if OK, false if canceled
*/
public
boolean
selectLMAParameters
(){
// int numSeries=fittingStrategy.getNumSeries();
boolean
resetCorrections
=
false
;
GenericDialog
gd
=
new
GenericDialog
(
"Levenberg-Marquardt algorithm parameters for cameras distortions/locations"
);
gd
.
addCheckbox
(
"Debug df/dX0, df/dY0"
,
false
);
gd
.
addCheckbox
(
"Keep current correction parameters"
,
this
.
keepCorrectionParameters
);
// gd.addNumericField("Iteration number to start (0.."+(numSeries-1)+")", this.seriesNumber, 0);
gd
.
addNumericField
(
"Initial LMA Lambda "
,
this
.
lambda
,
5
);
gd
.
addNumericField
(
"Multiply lambda on success"
,
this
.
lambdaStepDown
,
5
);
...
...
@@ -1192,6 +1188,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
gd
.
addCheckbox
(
"Show disabled parameters"
,
this
.
showDisabledParams
);
gd
.
addCheckbox
(
"Show per-sample correction parameters"
,
this
.
showCorrectionParams
);
gd
.
addCheckbox
(
"Reset all per-sample corrections to zero"
,
resetCorrections
);
// gd.addCheckbox("Show debug images before correction",this.showThisImages);
// gd.addCheckbox("Show debug images after correction", this.showNextImages);
...
...
@@ -1200,6 +1197,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
gd
.
showDialog
();
if
(
gd
.
wasCanceled
())
return
false
;
this
.
debugDerivativesFxDxDy
=
gd
.
getNextBoolean
();
this
.
keepCorrectionParameters
=
gd
.
getNextBoolean
();
// this.seriesNumber= (int) gd.getNextNumber();
this
.
lambda
=
gd
.
getNextNumber
();
this
.
lambdaStepDown
=
gd
.
getNextNumber
();
...
...
@@ -1218,10 +1216,12 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
// this.showNextImages= gd.getNextBoolean();
// this.threadsMax= (int) gd.getNextNumber();
// this.threadedLMA= gd.getNextBoolean();
resetCorrections
=
gd
.
getNextBoolean
();
if
(
resetCorrections
)
fieldFitting
.
resetSampleCorr
();
return
true
;
}
}
public
void
listParameters
(
String
title
,
String
path
){
public
void
listParameters
(
String
title
,
String
path
){
String
header
=
"State\tDescription\tValue\tUnits"
;
StringBuffer
sb
=
new
StringBuffer
();
for
(
String
s:
this
.
fieldFitting
.
getParameterValueStrings
(
true
,
true
)){
//this.showDisabledParams)){ //
...
...
@@ -1234,9 +1234,9 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
else
{
new
TextWindow
(
title
,
header
,
sb
.
toString
(),
800
,
1000
);
}
}
}
public
void
listData
(
String
title
,
String
path
){
public
void
listData
(
String
title
,
String
path
){
if
((
showSamples
==
null
)
||
(
showSamples
.
length
!=
sampleCoord
.
length
*
sampleCoord
[
0
].
length
)){
showSamples
=
new
boolean
[
sampleCoord
.
length
*
sampleCoord
[
0
].
length
];
for
(
int
i
=
0
;
i
<
showSamples
.
length
;
i
++)
showSamples
[
i
]=
false
;
...
...
@@ -1292,9 +1292,9 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
localShowSamples
,
showIgnoredData
,
showRad
);
}
}
public
void
listData
(
String
title
,
public
void
listData
(
String
title
,
String
path
,
boolean
showMotors
,
boolean
[]
showMeasCalc
,
...
...
@@ -1512,12 +1512,12 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
new
TextWindow
(
title
,
header
,
sb
.
toString
(),
800
,
1000
);
}
}
public
void
showCurvCorr
(){
}
public
void
showCurvCorr
(){
fieldFitting
.
showCurvCorr
(
"curv_corr"
);
}
}
public
void
listCombinedResults
(){
public
void
listCombinedResults
(){
// private boolean rslt_show_z_axial=true;
// private boolean rslt_show_z_individual=true;
// private boolean rslt_show_f_axial=true;
...
...
@@ -1583,9 +1583,9 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
center_z
[
1
]+
rslt_scan_above
,
//double scan_above,
rslt_scan_step
,
//double scan_step,
rslt_mtf50_mode
);
//boolean freq_mode)
}
}
public
void
listCombinedResults
(
public
void
listCombinedResults
(
String
title
,
String
path
,
boolean
[]
show_chn
,
...
...
@@ -1676,10 +1676,10 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
else
{
new
TextWindow
(
title
,
header
,
sb
.
toString
(),
800
,
1000
);
}
}
}
public
void
listScanQB
(){
public
void
listScanQB
(){
double
[]
center_z
=
fieldFitting
.
getZCenters
();
double
best_qb_axial
=
fieldFitting
.
getBestQualB
(
...
...
@@ -1725,9 +1725,9 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
this
.
k_blue
,
//double kb,
this
.
qb_use_corrected
,
//boolean corrected,
this
.
qb_invert
);
//boolean freq_mode)
}
}
public
void
listScanQB
(
public
void
listScanQB
(
String
title
,
String
path
,
double
min_z
,
// absolute
...
...
@@ -1766,10 +1766,10 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
else
{
new
TextWindow
(
title
,
header
,
sb
.
toString
(),
800
,
1000
);
}
}
}
public
boolean
dialogLMAStep
(
boolean
[]
state
){
public
boolean
dialogLMAStep
(
boolean
[]
state
){
String
[]
states
={
"Worse, increase lambda"
,
"Better, decrease lambda"
,
...
...
@@ -1836,11 +1836,11 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
// this.showNextImages= gd.getNextBoolean();
this
.
saveSeries
=
true
;
return
gd
.
wasOKed
();
}
}
public
boolean
LevenbergMarquardt
(
boolean
openDialog
,
int
debugLevel
){
public
boolean
LevenbergMarquardt
(
boolean
openDialog
,
int
debugLevel
){
double
savedLambda
=
this
.
lambda
;
this
.
debugLevel
=
debugLevel
;
if
(
openDialog
&&
!
selectLMAParameters
())
return
false
;
...
...
@@ -1940,7 +1940,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
// if (wasLastSeries) break;
// } // while (this.fittingStrategy.isSeriesValid(this.seriesNumber)){ // TODO: Add "stop" tag to series
if
(
debugLevel
>
0
)
System
.
out
.
println
(
"LevenbergMarquardt(): all done"
+
", step="
+
this
.
iterationStepNumber
+
", step="
+
this
.
iterationStepNumber
+
", RMS="
+
this
.
currentRMS
+
" ("
+
this
.
firstRMS
+
") "
+
" at "
+
IJ
.
d2s
(
0.000000001
*(
System
.
nanoTime
()-
this
.
startTime
),
3
));
...
...
@@ -1948,7 +1948,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
commitParameterVector
(
this
.
savedVector
);
return
true
;
// all series done
}
}
...
...
@@ -2116,6 +2116,8 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
hConfig
.
setThrowExceptionOnMissing
(
false
);
// default value, will return null on missing
comment
=
hConfig
.
getString
(
"comment"
,
"no comments"
);
// if ((comment.length()>10) && comment.substring(0,9).equals("<![CDATA[")) comment=comment.substring(9,comment.length()-3);
serialNumber
=
hConfig
.
getString
(
"serialNumber"
,
"???"
);
lensSerial
=
hConfig
.
getString
(
"lensSerial"
,
"???"
);
pX0_distortions
=
Double
.
parseDouble
(
hConfig
.
getString
(
"lens_center_x"
,
"0.0"
));
...
...
@@ -2154,7 +2156,11 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
String
path
){
// x,y,r
XMLConfiguration
hConfig
=
new
XMLConfiguration
();
hConfig
.
setRootElementName
(
"focusingHistory"
);
if
(
comment
!=
null
)
hConfig
.
addProperty
(
"comment"
,
comment
);
if
(
comment
!=
null
){
String
comment_esc
=
comment
.
replace
(
","
,
"\\,"
);
// hConfig.addProperty("comment","<![CDATA["+comment_esc+ "]]>");
hConfig
.
addProperty
(
"comment"
,
comment_esc
);
}
if
(
serialNumber
!=
null
)
hConfig
.
addProperty
(
"serialNumber"
,
serialNumber
);
if
(
lensSerial
!=
null
)
hConfig
.
addProperty
(
"lensSerial"
,
lensSerial
);
hConfig
.
addProperty
(
"lens_center_x"
,
pX0_distortions
);
// distortions center, not aberrations!
...
...
@@ -2206,7 +2212,10 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
private
double
[][]
sampleCorrPullZero
=
new
double
[
6
][];
// 1.0 - only difference from neighbors matters, 0.0 - only difference from 0
private
double
[]
sampleCorrRadius
=
null
;
private
double
[][][][]
sampleCorrCrossWeights
=
new
double
[
6
][][][];
private
double
[]
sampleCorr
=
null
;
private
double
[]
sampleCorrVector
=
null
;
// currently adjusted per-sample parameters
private
double
[][][]
correctionParameters
=
new
double
[
6
][][];
// all
public
int
numberOfLocations
=
0
;
private
int
[][]
sampleCorrChnParIndex
=
null
;
private
boolean
[]
dflt_sampleCorrSelect
=
{
false
,
false
,
false
,
false
};
private
double
[]
dflt_sampleCorrCost
=
{
0.2
,
50.0
,
1.0
,
1.0
};
...
...
@@ -2217,6 +2226,46 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
"Green, sagittal"
,
"Green, tangential"
,
"Blue, sagittal"
,
"Blue, tangential"
};
public
void
setProperties
(
String
prefix
,
Properties
properties
){
properties
.
setProperty
(
prefix
+
"centerSelect_X"
,
centerSelect
[
0
]+
""
);
properties
.
setProperty
(
prefix
+
"centerSelect_Y"
,
centerSelect
[
1
]+
""
);
mechanicalFocusingModel
.
setProperties
(
prefix
+
"mechanicalFocusingModel."
,
properties
);
for
(
int
i
=
0
;
i
<
curvatureModel
.
length
;
i
++){
if
(
curvatureModel
[
i
]!=
null
)
curvatureModel
[
i
].
setProperties
(
prefix
+
"curvatureModel_"
+
i
+
"."
,
properties
);
}
if
(
channelSelect
!=
null
)
for
(
int
i
=
0
;
i
<
channelSelect
.
length
;
i
++){
properties
.
setProperty
(
prefix
+
"channelSelect_"
+
i
,
channelSelect
[
i
]+
""
);
}
if
(
mechanicalSelect
!=
null
)
for
(
int
i
=
0
;
i
<
mechanicalSelect
.
length
;
i
++){
properties
.
setProperty
(
prefix
+
"mechanicalSelect_"
+
i
,
mechanicalSelect
[
i
]+
""
);
}
for
(
int
chn
=
0
;
chn
<
curvatureSelect
.
length
;
chn
++)
if
(
curvatureSelect
[
chn
]!=
null
)
for
(
int
i
=
0
;
i
<
curvatureSelect
[
chn
].
length
;
i
++){
properties
.
setProperty
(
prefix
+
"curvatureSelect_"
+
chn
+
"_"
+
i
,
curvatureSelect
[
chn
][
i
]+
""
);
}
for
(
int
chn
=
0
;
chn
<
sampleCorrSelect
.
length
;
chn
++)
if
(
sampleCorrSelect
[
chn
]!=
null
)
for
(
int
i
=
0
;
i
<
sampleCorrSelect
[
chn
].
length
;
i
++){
properties
.
setProperty
(
prefix
+
"sampleCorrSelect_"
+
chn
+
"_"
+
i
,
sampleCorrSelect
[
chn
][
i
]+
""
);
}
for
(
int
chn
=
0
;
chn
<
sampleCorrCost
.
length
;
chn
++)
if
(
sampleCorrCost
[
chn
]!=
null
)
for
(
int
i
=
0
;
i
<
sampleCorrCost
[
chn
].
length
;
i
++){
properties
.
setProperty
(
prefix
+
"sampleCorrCost_"
+
chn
+
"_"
+
i
,
sampleCorrCost
[
chn
][
i
]+
""
);
}
for
(
int
chn
=
0
;
chn
<
sampleCorrSigma
.
length
;
chn
++)
if
(
sampleCorrSigma
[
chn
]!=
null
)
for
(
int
i
=
0
;
i
<
sampleCorrSigma
[
chn
].
length
;
i
++){
properties
.
setProperty
(
prefix
+
"sampleCorrSigma_"
+
chn
+
"_"
+
i
,
sampleCorrSigma
[
chn
][
i
]+
""
);
}
for
(
int
chn
=
0
;
chn
<
sampleCorrPullZero
.
length
;
chn
++)
if
(
sampleCorrPullZero
[
chn
]!=
null
)
for
(
int
i
=
0
;
i
<
sampleCorrPullZero
[
chn
].
length
;
i
++){
properties
.
setProperty
(
prefix
+
"sampleCorrPullZero_"
+
chn
+
"_"
+
i
,
sampleCorrPullZero
[
chn
][
i
]+
""
);
}
/*
private double [][] sampleCorrCost= new double[6][]; // equivalent cost of one unit of parameter value (in result units, um)
private double [][] sampleCorrSigma= new double[6][]; // sigma (in mm) for neighbors influence
private double [][] sampleCorrPullZero=new double[6][]; // 1.0 - only difference from neighbors matters, 0.0 - only difference from 0
*/
}
public
double
[]
getSampleRadiuses
(){
// distance from the current center to each each sample
return
sampleCorrRadius
;
}
...
...
@@ -2224,24 +2273,18 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
int
width
=
getSampleWidth
();
int
numSamples
=
getNumSamples
();
String
[]
chnShortNames
={
"RS"
,
"RT"
,
"GS"
,
"GT"
,
"BRS"
,
"BT"
};
if
((
sampleCorr
==
null
)
||
(
sampleCorr
.
length
==
0
))
{
String
msg
=
"No correction parameters are enabled"
;
IJ
.
showMessage
(
msg
);
System
.
out
.
println
(
msg
);
return
;
}
int
numCorrPar
=
0
;
int
maxNumPars
=
0
;
for
(
int
chn
=
0
;
chn
<
sampleCorrChnParIndex
.
length
;
chn
++)
if
(
sampleCorrChnParIndex
[
chn
]!=
null
)
for
(
int
np
=
0
;
np
<
sampleCorrChnParIndex
[
chn
].
length
;
np
++)
if
(
sampleCorrChnParIndex
[
chn
][
np
]>=
0
)
{
for
(
int
chn
=
0
;
chn
<
correctionParameters
.
length
;
chn
++)
if
(
correctionParameters
[
chn
]!=
null
)
for
(
int
np
=
0
;
np
<
correctionParameters
[
chn
].
length
;
np
++)
if
(
correctionParameters
[
chn
][
np
]!=
null
)
{
numCorrPar
++;
if
(
np
>
maxNumPars
)
maxNumPars
=
np
;
}
maxNumPars
++;
if
(
numCorrPar
==
0
){
String
msg
=
"Still n
o correction parameters are enabled"
;
String
msg
=
"N
o correction parameters are enabled"
;
IJ
.
showMessage
(
msg
);
System
.
out
.
println
(
msg
);
return
;
...
...
@@ -2249,25 +2292,12 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
double
[][]
pixels
=
new
double
[
numCorrPar
][
numSamples
];
String
[]
titles
=
new
String
[
numCorrPar
];
int
index
=
0
;
/*
for (int chn=0;chn<sampleCorrChnParIndex.length;chn++)
if (sampleCorrChnParIndex[chn]!=null)
for (int np=0;np<sampleCorrChnParIndex[chn].length;np++)
if (sampleCorrChnParIndex[chn][np]>=0) {
titles[index]=chnShortNames[chn]+"-"+curvatureModel[chn].getZName(np);
for (int nSample=0;nSample<numSamples;nSample++) {
pixels[index][nSample]=sampleCorr[sampleCorrChnParIndex[chn][np]+nSample];
}
index++;
}
*/
for
(
int
np
=
0
;
np
<
maxNumPars
;
np
++)
for
(
int
chn
=
0
;
chn
<
sampleCorrChnParIndex
.
length
;
chn
++)
if
((
sampleCorrChnParIndex
[
chn
]!=
null
)
&&
(
sampleCorrChnParIndex
[
chn
].
length
>
np
)
&&
(
sampleCorrChnParIndex
[
chn
][
np
]>=
0
))
{
for
(
int
chn
=
0
;
chn
<
correctionParameters
.
length
;
chn
++)
if
((
correctionParameters
[
chn
]!=
null
)
&&
(
correctionParameters
[
chn
].
length
>
np
)
&&
(
correctionParameters
[
chn
][
np
]!=
null
))
{
titles
[
index
]=
chnShortNames
[
chn
]+
"-"
+
curvatureModel
[
chn
].
getZName
(
np
);
for
(
int
nSample
=
0
;
nSample
<
numSamples
;
nSample
++)
{
pixels
[
index
][
nSample
]=
sampleCorr
[
sampleCorrChnParIndex
[
chn
][
np
]+
nSample
];
pixels
[
index
][
nSample
]=
correctionParameters
[
chn
][
np
][
nSample
];
}
index
++;
}
...
...
@@ -2489,7 +2519,10 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
int
fromChn
=
individualChannels
?
0
:
firstChn
;
int
toChn
=
individualChannels
?
channelSelect
.
length
:(
firstChn
+
1
);
boolean
resetCorrections
=
false
;
GenericDialog
gd
=
new
GenericDialog
(
title
);
gd
.
addCheckbox
(
"Reset all per-sample corrections to zero"
,
resetCorrections
);
int
numParsZR
[]
=
getNumCurvars
();
// [0] - Z, [1] - r,
for
(
int
nChn
=
fromChn
;
nChn
<
toChn
;
nChn
++)
if
(
channelSelect
[
nChn
]){
...
...
@@ -2507,6 +2540,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
gd
.
showDialog
();
if
(
gd
.
wasCanceled
())
return
false
;
if
(
gd
.
wasOKed
())
{
// selected non-default "Apply"
resetCorrections
=
gd
.
getNextBoolean
();
for
(
int
nChn
=
fromChn
;
nChn
<
toChn
;
nChn
++)
if
(
channelSelect
[
nChn
]){
for
(
int
i
=
0
;
i
<
sampleCorrSelect
[
nChn
].
length
;
i
++)
if
(
disabledPars
||
curvatureSelect
[
nChn
][
i
*
numParsZR
[
1
]]){
sampleCorrSelect
[
nChn
][
i
]=
gd
.
getNextBoolean
();
...
...
@@ -2528,27 +2562,86 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
}
}
if
(
resetCorrections
)
resetSampleCorr
();
return
true
;
}
// once per data set
public
void
resetSampleCorr
(){
for
(
int
i
=
0
;
i
<
correctionParameters
.
length
;
i
++)
correctionParameters
[
i
]=
null
;
}
public
double
[]
getCorrVector
(){
//numberOfLocations
int
numPars
=
0
;
for
(
int
nChn
=
0
;
nChn
<
sampleCorrChnParIndex
.
length
;
nChn
++)
{
if
(
sampleCorrChnParIndex
[
nChn
]!=
null
){
for
(
int
nPar
=
0
;
nPar
<
sampleCorrChnParIndex
[
nChn
].
length
;
nPar
++)
{
if
(
sampleCorrChnParIndex
[
nChn
][
nPar
]>=
0
){
numPars
+=
numberOfLocations
;
}
}
}
}
sampleCorrVector
=
new
double
[
numPars
];
for
(
int
nChn
=
0
;
nChn
<
sampleCorrChnParIndex
.
length
;
nChn
++)
{
if
(
sampleCorrChnParIndex
[
nChn
]!=
null
){
for
(
int
nPar
=
0
;
nPar
<
sampleCorrChnParIndex
[
nChn
].
length
;
nPar
++)
{
if
(
sampleCorrChnParIndex
[
nChn
][
nPar
]>=
0
){
for
(
int
i
=
0
;
i
<
numberOfLocations
;
i
++){
if
((
correctionParameters
[
nChn
]!=
null
)
&&
(
correctionParameters
[
nChn
][
nPar
]!=
null
))
sampleCorrVector
[
sampleCorrChnParIndex
[
nChn
][
nPar
]+
i
]=
correctionParameters
[
nChn
][
nPar
][
i
];
else
sampleCorrVector
[
sampleCorrChnParIndex
[
nChn
][
nPar
]+
i
]=
0.0
;
}
}
}
}
}
return
sampleCorrVector
;
}
public
void
commitCorrVector
(){
commitCorrVector
(
sampleCorrVector
);
}
public
void
commitCorrVector
(
double
[]
vector
){
for
(
int
nChn
=
0
;
nChn
<
sampleCorrChnParIndex
.
length
;
nChn
++)
{
if
(
sampleCorrChnParIndex
[
nChn
]!=
null
){
for
(
int
nPar
=
0
;
nPar
<
sampleCorrChnParIndex
[
nChn
].
length
;
nPar
++)
{
if
(
sampleCorrChnParIndex
[
nChn
][
nPar
]>=
0
){
if
(
correctionParameters
[
nChn
]==
null
)
correctionParameters
[
nChn
]=
new
double
[
sampleCorrChnParIndex
[
nChn
].
length
][];
if
(
correctionParameters
[
nChn
][
nPar
]==
null
)
correctionParameters
[
nChn
][
nPar
]=
new
double
[
numberOfLocations
];
for
(
int
i
=
0
;
i
<
numberOfLocations
;
i
++){
correctionParameters
[
nChn
][
nPar
][
i
]=
vector
[
sampleCorrChnParIndex
[
nChn
][
nPar
]+
i
];
}
}
}
}
}
}
/**
* create matrix of weights of the other parameters influence
* @param sampleCoordinates [sample number]{x,y} - flattened array of sample coordinates
* Run in the beginning of fitting series (zeroes the values)
*/
public
void
initSampleCorr
(
double
[][]
sampleCoordinates
){
int
numSamples
=
sampleCoordinates
.
length
;
// once per fitting series (or parameter change
public
void
initSampleCorrVector
(
double
[][]
sampleCoordinates
){
numberOfLocations
=
sampleCoordinates
.
length
;
// int numSamples=sampleCoordinates.length;
for
(
int
nChn
=
0
;
nChn
<
sampleCorrSelect
.
length
;
nChn
++)
{
if
(
channelSelect
[
nChn
]){
int
numPars
=
sampleCorrSelect
[
nChn
].
length
;
sampleCorrCrossWeights
[
nChn
]=
new
double
[
numPars
][][];
for
(
int
nPar
=
0
;
nPar
<
numPars
;
nPar
++){
if
(
sampleCorrSelect
[
nChn
][
nPar
]){
sampleCorrCrossWeights
[
nChn
][
nPar
]=
new
double
[
numSamples
][
numSample
s
];
sampleCorrCrossWeights
[
nChn
][
nPar
]=
new
double
[
numberOfLocations
][
numberOfLocation
s
];
double
k
=-
getPixelMM
()*
getPixelMM
()/(
sampleCorrSigma
[
nChn
][
nPar
]*
sampleCorrSigma
[
nChn
][
nPar
]);
for
(
int
i
=
0
;
i
<
numSample
s
;
i
++){
for
(
int
i
=
0
;
i
<
numberOfLocation
s
;
i
++){
double
sw
=
0.0
;
for
(
int
j
=
0
;
j
<
numSample
s
;
j
++)
{
for
(
int
j
=
0
;
j
<
numberOfLocation
s
;
j
++)
{
if
(
i
!=
j
){
double
dx
=
sampleCoordinates
[
i
][
0
]-
sampleCoordinates
[
i
][
0
];
double
dy
=
sampleCoordinates
[
i
][
1
]-
sampleCoordinates
[
i
][
1
];
...
...
@@ -2559,7 +2652,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
if
((
sampleCorrPullZero
[
nChn
][
nPar
]==
0
)
||(
sampleCorrCost
[
nChn
][
nPar
]==
0
))
sw
=
0.0
;
else
if
(
sw
!=
0.0
)
sw
=-
sampleCorrCost
[
nChn
][
nPar
]*
sampleCorrPullZero
[
nChn
][
nPar
]/
sw
;
for
(
int
j
=
0
;
j
<
numSample
s
;
j
++)
{
for
(
int
j
=
0
;
j
<
numberOfLocation
s
;
j
++)
{
if
(
i
!=
j
){
sampleCorrCrossWeights
[
nChn
][
nPar
][
i
][
j
]*=
sw
;
}
else
{
...
...
@@ -2593,11 +2686,12 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
}
// currently all correction parameters are initialized as zeros.
sampleCorr
=
new
double
[
numPars
];
for
(
int
i
=
0
;
i
<
numPars
;
i
++)
sampleCorr
[
i
]=
0.0
;
sampleCorrRadius
=
new
double
[
numSamples
];
getCorrVector
();
// sampleCorrVector=new double [numPars];
// for (int i=0;i<numPars;i++)sampleCorrVector[i]=0.0;
sampleCorrRadius
=
new
double
[
numberOfLocations
];
//pXY
for
(
int
i
=
0
;
i
<
numSample
s
;
i
++){
for
(
int
i
=
0
;
i
<
numberOfLocation
s
;
i
++){
double
dx
=
sampleCoordinates
[
i
][
0
]-
pXY
[
0
];
double
dy
=
sampleCoordinates
[
i
][
1
]-
pXY
[
0
];
sampleCorrRadius
[
i
]=
getPixelMM
()*
Math
.
sqrt
(
dx
*
dx
+
dy
*
dy
);
...
...
@@ -2605,12 +2699,20 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
public
double
[]
getCorrPar
(
int
chn
,
int
sampleIndex
){
if
((
sampleCorrChnParIndex
==
null
)
||
(
sampleCorrChnParIndex
[
chn
]==
null
))
return
null
;
/*
if ((sampleCorrChnParIndex==null) || (sampleCorrChnParIndex[chn]==null)) return null;
double [] corr =new double [sampleCorrChnParIndex[chn].length];
for (int i=0;i<corr.length;i++){
if (sampleCorrChnParIndex[chn][i]<0) corr[i]=0.0;
else
corr
[
i
]=
sampleCor
r
[
sampleCorrChnParIndex
[
chn
][
i
]+
sampleIndex
];
else corr[i]=sampleCorrVecto
r[sampleCorrChnParIndex[chn][i]+sampleIndex];
}
*/
if
(
correctionParameters
[
chn
]==
null
)
return
null
;
double
[]
corr
=
new
double
[
correctionParameters
[
chn
].
length
];
for
(
int
i
=
0
;
i
<
corr
.
length
;
i
++){
if
(
correctionParameters
[
chn
][
i
]
!=
null
)
corr
[
i
]=
correctionParameters
[
chn
][
i
][
sampleIndex
];
else
corr
[
i
]=
0.0
;
}
return
corr
;
}
...
...
@@ -2751,7 +2853,8 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
disabledCorrectionPars
))
return
false
;
// initSampleCorr(flattenSampleCoord());
}
initSampleCorr
(
flattenSampleCoord
());
// run always regardless of configured or not (to create zero-length array of corr)
// will modify
initSampleCorrVector
(
flattenSampleCoord
());
// run always regardless of configured or not (to create zero-length array of corr)
return
true
;
}
ArrayList
<
String
>
getParameterValueStrings
(
boolean
showDisabled
,
boolean
showCorrection
){
...
...
@@ -2804,12 +2907,12 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
if
(
showCorrection
&&
(
getNumberOfCorrParameters
()>
0
)){
parList
.
add
(
"\t ===== Per-sample correction parameters =====\t\t"
);
for
(
int
n
=
0
;
n
<
sampleCorrChnParIndex
.
length
;
n
++)
if
(
sampleCorrChnParIndex
[
n
]!=
null
){
for
(
int
np
=
0
;
np
<
sampleCorrChnParIndex
[
n
].
length
;
np
++)
if
(
sampleCorrChnParIndex
[
n
][
np
]>=
0
){
int
numSamples
=
sampleCorrCrossWeights
[
n
][
np
].
length
;
for
(
int
n
=
0
;
n
<
correctionParameters
.
length
;
n
++)
if
(
correctionParameters
[
n
]!=
null
){
for
(
int
np
=
0
;
np
<
correctionParameters
[
n
].
length
;
np
++)
if
(
correctionParameters
[
n
][
np
]!=
null
){
//
int numSamples=sampleCorrCrossWeights[n][np].length;
parList
.
add
(
"\t ----- correction parameters for \""
+
getDescription
(
n
)+
" "
+
curvatureModel
[
n
].
getZDescription
(
np
)+
"\" -----\t\t"
);
for
(
int
i
=
0
;
i
<
numSample
s
;
i
++){
parList
.
add
(
i
+
"\t"
+
curvatureModel
[
n
].
getZDescription
(
np
)+
":\t"
+
sampleCorr
[
sampleCorrChnParIndex
[
n
][
np
]+
i
]+
"\t"
);
for
(
int
i
=
0
;
i
<
numberOfLocation
s
;
i
++){
parList
.
add
(
i
+
"\t"
+
curvatureModel
[
n
].
getZDescription
(
np
)+
":\t"
+
correctionParameters
[
n
][
np
][
i
]+
"\t"
);
}
}
}
...
...
@@ -2892,8 +2995,8 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
return
true
;
}
public
int
getNumberOfCorrParameters
(){
return
(
sampleCorr
!=
null
)?
sampleCor
r
.
length
:
0
;
public
int
getNumberOfCorrParameters
(){
// selected for fitting
return
(
sampleCorrVector
!=
null
)?
sampleCorrVecto
r
.
length
:
0
;
}
public
int
getNumberOfParameters
(
boolean
sagittalMaster
){
...
...
@@ -2954,7 +3057,7 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
}
int
nCorrPars
=
getNumberOfCorrParameters
();
for
(
int
i
=
0
;
i
<
nCorrPars
;
i
++)
pars
[
np
++]=
sampleCor
r
[
i
];
for
(
int
i
=
0
;
i
<
nCorrPars
;
i
++)
pars
[
np
++]=
sampleCorrVecto
r
[
i
];
return
pars
;
}
...
...
@@ -2982,7 +3085,8 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
}
// copy correction parameters
int
nCorrPars
=
getNumberOfCorrParameters
();
for
(
int
i
=
0
;
i
<
nCorrPars
;
i
++)
sampleCorr
[
i
]=
pars
[
np
++];
for
(
int
i
=
0
;
i
<
nCorrPars
;
i
++)
sampleCorrVector
[
i
]=
pars
[
np
++];
commitCorrVector
();
// copy center parameters to dependent
// copy if master is selected, regardless of is dependent selected or not
...
...
@@ -3197,6 +3301,22 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
public
MechanicalFocusingModel
(){
// add arguments?
initDefaults
();
}
public
void
setProperties
(
String
prefix
,
Properties
properties
){
if
(
paramValues
==
null
)
{
for
(
int
i
=
0
;
i
<
paramValues
.
length
;
i
++){
properties
.
setProperty
(
prefix
+
descriptions
[
i
][
0
],
paramValues
[
i
]+
""
);
}
}
}
public
void
getProperties
(
String
prefix
,
Properties
properties
){
if
((
paramValues
==
null
)
||
(
paramValues
.
length
!=
descriptions
.
length
))
initDefaults
();
for
(
int
i
=
0
;
i
<
paramValues
.
length
;
i
++){
if
(
properties
.
getProperty
(
prefix
+
descriptions
[
i
][
0
])!=
null
)
paramValues
[
i
]=
Double
.
parseDouble
(
properties
.
getProperty
(
prefix
+
descriptions
[
i
][
0
]));
}
}
public
void
initDefaults
(){
paramValues
=
new
double
[
descriptions
.
length
];
for
(
int
i
=
0
;
i
<
descriptions
.
length
;
i
++)
paramValues
[
i
]=
Double
.
parseDouble
(
descriptions
[
i
][
3
]);
...
...
@@ -3456,6 +3576,46 @@ d_s2/d_x0= 2*delta_x*delta_y^2/r2^2
this
.
modelParams
=
new
double
[
distanceParametersNumber
][
radialParametersNumber
];
setDefaults
();
}
public
void
setProperties
(
String
prefix
,
Properties
properties
){
properties
.
setProperty
(
prefix
+
"distanceParametersNumber"
,
modelParams
.
length
+
""
);
properties
.
setProperty
(
prefix
+
"radialParametersNumber"
,
modelParams
[
0
].
length
+
""
);
properties
.
setProperty
(
prefix
+
"pX0"
,
pX0
+
""
);
properties
.
setProperty
(
prefix
+
"pY0"
,
pY0
+
""
);
for
(
int
i
=
0
;
i
<
modelParams
.
length
;
i
++)
for
(
int
j
=
0
;
j
<
modelParams
[
i
].
length
;
j
++){
properties
.
setProperty
(
prefix
+
"modelParams"
+
i
+
"_"
+
j
,
modelParams
[
i
][
j
]+
""
);
}
}
public
void
getProperties
(
String
prefix
,
Properties
properties
){
int
numZPars
,
numRPars
;
if
(
properties
.
getProperty
(
prefix
+
"pX0"
)!=
null
)
pX0
=
Double
.
parseDouble
(
properties
.
getProperty
(
prefix
+
"pX0"
));
if
(
properties
.
getProperty
(
prefix
+
"pY0"
)!=
null
)
pX0
=
Double
.
parseDouble
(
properties
.
getProperty
(
prefix
+
"pY0"
));
if
(
modelParams
!=
null
){
numZPars
=
modelParams
.
length
;
numRPars
=
modelParams
[
0
].
length
;
}
else
{
numZPars
=
dflt_distanceParametersNumber
;
numRPars
=
dflt_radialParametersNumber
;
}
if
(
properties
.
getProperty
(
prefix
+
"distanceParametersNumber"
)!=
null
)
numZPars
=
Integer
.
parseInt
(
properties
.
getProperty
(
prefix
+
"distanceParametersNumber"
));
if
(
properties
.
getProperty
(
prefix
+
"radialParametersNumber"
)!=
null
)
numRPars
=
Integer
.
parseInt
(
properties
.
getProperty
(
prefix
+
"radialParametersNumber"
));
if
((
modelParams
!=
null
)
||
(
modelParams
.
length
!=
numZPars
)||
(
modelParams
[
0
].
length
!=
numRPars
)){
modelParams
=
new
double
[
numZPars
][
numRPars
];
setDefaults
();
}
for
(
int
i
=
0
;
i
<
modelParams
.
length
;
i
++)
for
(
int
j
=
0
;
j
<
modelParams
[
i
].
length
;
j
++){
if
(
properties
.
getProperty
(
prefix
+
"modelParams"
+
i
+
"_"
+
j
)!=
null
)
{
modelParams
[
i
][
j
]=
Double
.
parseDouble
(
properties
.
getProperty
(
prefix
+
"modelParams"
+
i
+
"_"
+
j
));
}
}
}
public
void
setCenter
(
double
pX
,
double
pY
){
pX0
=
pX
;
pY0
=
pY
;
...
...
@@ -3772,3 +3932,6 @@ ar1 - ar[4]
}
}
}
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment