inaccurate result in function polyfit (Bug #2887)

Added by Zhoucheng Li almost 12 years ago. Updated over 11 years ago.

Status:Done Start date:2013-03-14
Priority:Normal Due date:
Assignee:Zhoucheng Li % Done:

0%
Category:legacy, contrib
Target version:3.0
Affected version:branch 'master' (3.0-dev) Operating System:Any
Difficulty: HW Platform:Any
Pull request:https://github.com/Itseez/opencv/pull/674

Description

The default datatype in polyfit is CV_32F and the method is Least square method.
In some circumstances where input data on x-axis are uniformly large enough(e.g. 24-47, 24 input numbers) while data on y-axis are relatively close and large too, the fitting curve could be tremendously inaccurate.

Example(pseudocode):

x1 = [0;1;2;3;4;5;6;7;8;9;10;11;12;13;14;15;16;17;18;19;20;21;22;23];
y1 = [488.3692016601563
488.4369506835938
488.4021911621094
488.3995666503906
488.3041076660156
488.3595581054688
488.2347412109375
488.2492065429688
488.0935974121094
488.1622314453125
488.1072082519531
488.2424011230469
488.1407165527344
488.2175903320313
488.2396850585938
488.27587890625];

x2 = [24;25;26;27;28;29;30;31;32;33;34;35;36;37;38;39;40;41;42;43;44;45;46;47];
y2 = [488.3522644042969
488.2658996582031
488.3119812011719
488.3011474609375
488.3711853027344
488.3617248535156
488.3223266601563
488.2750549316406
488.2916259765625
488.3085632324219
488.3853454589844
488.3346862792969
488.2177734375
488.3077392578125
488.29150390625
488.3375854492188];

Mat result1, result2;
polyfit(x1, y1, result1, 2);//fitting a second order curve
polyfit(x2, y2, result2, 2);//fitting a second order curve

Result(datatype:CV_32F, see picture "datatype : CV_32F"):

result1 = [488.4185180664063
-0.03082548454403877
0.001051958999596536];

result2 = [488.6778259277344
-0.02593995630741119
0.0003619703638833016];

I modified the code to use datatyoe CV_64F.
X = Mat::zeros(src_x.rows, order+1,CV_64F);
Mat srcy;
src_y.convertTo(srcy,CV_64F);
for(int i = 0; i <=order;i++){
src_x.convertTo(copy,CV_64F);
pow(copy,i,copy);
Mat M1 = X.col(i);
copy.col(0).copyTo(M1);}

Result(datatype:CV_64F, see picture "datatype : CV_64F" ):

result1 = [488.4181859764685
-0.03089823452888219
0.001054061977798282];

result2 = [488.6890140021986
-0.02416449552013233
0.0003906240721791021];

I think this is because for some input data on x-axis, the sum of x.^4 could be extremly large which may result in a extremly inaccurate reverse matrix.(In the example above, the sum of x.^4 is 46912204)

CV_32F.png - datatype : CV_32F (5.8 kB) Zhoucheng Li, 2013-03-14 03:38 pm

CV_64F.png - datatype : CV_64F (5.8 kB) Zhoucheng Li, 2013-03-14 03:38 pm

Associated revisions

Revision 6d4c4dcd
Added by Alexander Alekhin over 10 years ago

Merge pull request #2887 from ilya-lavrenov:ipp_morph_fix

History

#1
Updated by Zhoucheng Li almost 12 years ago

Another way to avoid this problem is to pre-process the data before using the function(datatype : CV_32F).
One can simply do the following procedures(The example is order 2)
  1. foreach element x2(n) and y2(n) in x2 and y2 respectively
    x2(n) -= average_x;
    y2(n) -= average_y;
  2. polyfit(x2, y2, result2, 2);
  3. Since the two curves have the same shape but are at different locations,we can derive the ralations between these 2 sets of curve parameters.
    Take "y = a*x.^2 + b*x + c" as an example.The new curve we've got is "y - average_y = a'*(x - average_x)^2 + b'*(x - average_x) + c'".
    From step2, we've got parameters a',b'and c'.
    Compare these 2 equations, we can derive that
    a = a';
    b = b' - 2*average_x*a';
    c = c' + average_y + a'*(average_x)^2 - b'*average_x;
  4. thus
    new_result2 = [488.6890572353125
    -0.024164495611330
    0.0003906240744981915];
    The answer is very close to the CV_64F version one.
    result2 = [488.6890140021986
    -0.02416449552013233
    0.0003906240721791021];

However, with a higher order(>=4),it's gonna be tough to derive the relations between parameters by paper work.

#2
Updated by Daniil Osokin almost 12 years ago

Hi, Zhoucheng! Could you contribute a bugfix (http://code.opencv.org/projects/opencv/wiki/How_to_contribute). It will speedup fixing in many times! Anyway, thank you for attention.

  • Assignee changed from Vadim Pisarevsky to Zhoucheng Li

#3
Updated by Zhoucheng Li almost 12 years ago

Daniil Osokin wrote:

Hi, Zhoucheng! Could you contribute a bugfix (http://code.opencv.org/projects/opencv/wiki/How_to_contribute). It will speedup fixing in many times! Anyway, thank you for attention.

My pleasure :D

#4
Updated by Daniil Osokin almost 12 years ago

Great, thanks!

  • Pull request set to https://github.com/Itseez/opencv/pull/674

#5
Updated by Andrey Kamaev almost 12 years ago

  • Target version changed from 2.4.5 to 3.0

#6
Updated by Dmitry Retinskiy over 11 years ago

  • Status changed from Open to Done
  • HW Platform set to Any
  • Operating System set to Any
  • Affected version changed from 2.4.4 to branch 'master' (3.0-dev)

Also available in: Atom PDF