OpenCV Q&A Forum - RSS feedhttp://answers.opencv.org/questions/OpenCV answersenCopyright <a href="http://www.opencv.org">OpenCV foundation</a>, 2012-2018.Wed, 16 Oct 2019 13:52:32 -0500Unclear how calibrateCamera estimates stdDeviations (perhaps wrong)http://answers.opencv.org/question/219805/unclear-how-calibratecamera-estimates-stddeviations-perhaps-wrong/Hi, I'm researching camera calibration and uncertainty propagation. I'm trying to understand how calibrateCameraExtended estimates stdDeviationsIntrinsics and stdDeviationsExtrinsics. [The docs](https://docs.opencv.org/ref/master/d9/d0c/group__calib3d.html#ga3207604e4b1a1758aa66acb6ed5aa65d) say very little.
So I go straight to the source code. Check [this lines of code](https://github.com/opencv/opencv/blob/00f36a261b09d6d6a89d41fb6268bceecd2e0dfd/modules/calib3d/src/calibration.cpp#L1798) where I understand the calculation is made.
It starts with `sigma2` (the "deviation of the noise") calculated as `norm(allErrors, NORM_L2SQR) / (total - nparams_nz);`
which is just the formula for the [unbiased estimator of the variance](https://en.wikipedia.org/wiki/Unbiased_estimation_of_standard_deviation). Ok so far.
And then it calculates each `s`-element of the vector of standard deviations `stdDevsM` by
stdDevsM.at<double>(s) = std::sqrt(JtJinv.at<double>(j,j) * sigma2);
Where `JtJinv` is the pseudo-inverse of the jacobian calculated a few lines above from `_JtJ` which in turn comes from the
LM solver invoked in [previous lines](https://github.com/opencv/opencv/blob/00f36a261b09d6d6a89d41fb6268bceecd2e0dfd/modules/calib3d/src/calibration.cpp#L1689).
**First question:** what exactly is `_JtJ`? I assume it must be the 1xN Jacobian of the projection error with respect to the parameters (there are N parameters). I've tried to trace the calculation of all the way to its origin, I got
[this far](https://github.com/opencv/opencv/blob/00f36a261b09d6d6a89d41fb6268bceecd2e0dfd/modules/calib3d/src/calibration.cpp#L809)
, but I'm not sure.
**Second question:** The moore-penrose of a 1xN matrix is a Nx1 matrix. So calling `JtJinv.at<double>(j,j)` with two indices `j` confuses me.
**Third question:** I couldn't make sense of the formula itself and I think It's wrong. The code has the comment
//see any papers about variance of the least squares estimator for
//detailed description of the variance estimation methods
But my understanding from simple [uncertainty propagation](https://en.wikipedia.org/wiki/Propagation_of_uncertainty#Linear_combinations) for the case of uncorrelated parameters is that the Jacobian vector `_JtJ`, the variance of the projection error `sigma2` and the vector of parameters standard deviations `stdDevsM` should follow (in matlab-like pseudocode):
sigma2 = dot(_JtJ.^2, stdDevsM.^2)
Which is ill-conditioned, the simplest solution would be
stdDevsM = sqrt(sigma2) .*_JtJ ./ norm(_JtJ.^2)
Is my reasoning correct? and, where does the calculation implemented in OpenCV come from?
Tue, 15 Oct 2019 20:26:26 -0500http://answers.opencv.org/question/219805/unclear-how-calibratecamera-estimates-stddeviations-perhaps-wrong/Answer by Eduardo for <p>Hi, I'm researching camera calibration and uncertainty propagation. I'm trying to understand how calibrateCameraExtended estimates stdDeviationsIntrinsics and stdDeviationsExtrinsics. <a href="https://docs.opencv.org/ref/master/d9/d0c/group__calib3d.html#ga3207604e4b1a1758aa66acb6ed5aa65d">The docs</a> say very little.</p>
<p>So I go straight to the source code. Check <a href="https://github.com/opencv/opencv/blob/00f36a261b09d6d6a89d41fb6268bceecd2e0dfd/modules/calib3d/src/calibration.cpp#L1798">this lines of code</a> where I understand the calculation is made.
It starts with <code>sigma2</code> (the "deviation of the noise") calculated as <code>norm(allErrors, NORM_L2SQR) / (total - nparams_nz);</code>
which is just the formula for the <a href="https://en.wikipedia.org/wiki/Unbiased_estimation_of_standard_deviation">unbiased estimator of the variance</a>. Ok so far.</p>
<p>And then it calculates each <code>s</code>-element of the vector of standard deviations <code>stdDevsM</code> by</p>
<pre><code>stdDevsM.at<double>(s) = std::sqrt(JtJinv.at<double>(j,j) * sigma2);
</code></pre>
<p>Where <code>JtJinv</code> is the pseudo-inverse of the jacobian calculated a few lines above from <code>_JtJ</code> which in turn comes from the
LM solver invoked in <a href="https://github.com/opencv/opencv/blob/00f36a261b09d6d6a89d41fb6268bceecd2e0dfd/modules/calib3d/src/calibration.cpp#L1689">previous lines</a>.</p>
<p><strong>First question:</strong> what exactly is <code>_JtJ</code>? I assume it must be the 1xN Jacobian of the projection error with respect to the parameters (there are N parameters). I've tried to trace the calculation of all the way to its origin, I got
<a href="https://github.com/opencv/opencv/blob/00f36a261b09d6d6a89d41fb6268bceecd2e0dfd/modules/calib3d/src/calibration.cpp#L809">this far</a>
, but I'm not sure.</p>
<p><strong>Second question:</strong> The moore-penrose of a 1xN matrix is a Nx1 matrix. So calling <code>JtJinv.at<double>(j,j)</code> with two indices <code>j</code> confuses me.</p>
<p><strong>Third question:</strong> I couldn't make sense of the formula itself and I think It's wrong. The code has the comment</p>
<pre><code>//see any papers about variance of the least squares estimator for
//detailed description of the variance estimation methods
</code></pre>
<p>But my understanding from simple <a href="https://en.wikipedia.org/wiki/Propagation_of_uncertainty#Linear_combinations">uncertainty propagation</a> for the case of uncorrelated parameters is that the Jacobian vector <code>_JtJ</code>, the variance of the projection error <code>sigma2</code> and the vector of parameters standard deviations <code>stdDevsM</code> should follow (in matlab-like pseudocode):</p>
<pre><code>sigma2 = dot(_JtJ.^2, stdDevsM.^2)
</code></pre>
<p>Which is ill-conditioned, the simplest solution would be</p>
<pre><code>stdDevsM = sqrt(sigma2) .*_JtJ ./ norm(_JtJ.^2)
</code></pre>
<p>Is my reasoning correct? and, where does the calculation implemented in OpenCV come from?</p>
http://answers.opencv.org/question/219805/unclear-how-calibratecamera-estimates-stddeviations-perhaps-wrong/?answer=219820#post-id-219820Quick answer:
1. `J` is the Jacobian that should be of size `MxN` with `M` the number of samples and `N` the number of estimated parameters; `JtJ` is simply `J^T x J` that should be of size `NxN`
2. See 1
3. I believe the [following](https://stat.ethz.ch/~geer/bsa199_o.pdf) applied to [non-linear least-squares](https://see.stanford.edu/materials/lsoeldsee263/07-ls-reg.pdf)/non-linear minimization should have the necessary informationWed, 16 Oct 2019 03:23:44 -0500http://answers.opencv.org/question/219805/unclear-how-calibratecamera-estimates-stddeviations-perhaps-wrong/?answer=219820#post-id-219820Comment by JamesBaxter for <p>Quick answer:</p>
<ol>
<li><code>J</code> is the Jacobian that should be of size <code>MxN</code> with <code>M</code> the number of samples and <code>N</code> the number of estimated parameters; <code>JtJ</code> is simply <code>J^T x J</code> that should be of size <code>NxN</code></li>
<li>See 1</li>
<li>I believe the <a href="https://stat.ethz.ch/~geer/bsa199_o.pdf">following</a> applied to <a href="https://see.stanford.edu/materials/lsoeldsee263/07-ls-reg.pdf">non-linear least-squares</a>/non-linear minimization should have the necessary information</li>
</ol>
http://answers.opencv.org/question/219805/unclear-how-calibratecamera-estimates-stddeviations-perhaps-wrong/?comment=219868#post-id-219868Perfecto, thank you! I spent a long time trying to trace the meaning of `JtJ` and couldn't figure it out. Also I was confusing the Jacobian wrt the optimisation function (projection error) with the Jacobian of the projection function.Wed, 16 Oct 2019 13:52:32 -0500http://answers.opencv.org/question/219805/unclear-how-calibratecamera-estimates-stddeviations-perhaps-wrong/?comment=219868#post-id-219868