curve fitting... tell me about it!!! :( [MATLAB]

[marcelino]

minsan ko nang naicharacterized the CM-R humidity sensor. pero ginamitan ko lang ng software. however, in the end kakailanganin ko pa din yung equation of how will i get it.

eto yung resistance curve ng humidity sensor.


tapos i am using 10-bits of adc to receive the resistance in terms of voltage... i just get the percentage of the humidity 100% as 1024. so i have a range of 0-1024[0-100%]. from 20-90% lang naman ang kailangan ko kasi yun lang ang response na valid for CM-R humidity sensor.

so eto na yun na generate kong curve fitted value (in terms of voltage (0-1024,0-5V)). ok naman in yung response, kinocompare ko sa isang electronic na humidity meter. BTW, ang humidty is actually a function of temperature.. kaya for the table below, may effect ang temperature sa humidity.

Code:

humidity (20 - 90%)
temperature (15 - 35 deg)

15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35
20 4 4 4 5 5 6 6 6 7 7 8 8 9 9 10 10 11 11 12 12 13
21 4 5 5 5 6 6 7 7 8 8 9 9 10 10 11 11 12 12 13 13 14
22 5 5 5 6 6 7 7 8 8 9 10 10 11 11 12 13 13 14 14 15 15
23 5 6 6 7 7 8 8 9 9 10 11 11 12 13 13 14 15 15 16 17 17
24 6 7 7 7 8 9 9 10 10 11 12 13 13 14 15 16 16 17 18 19 20
25 7 8 8 9 9 10 10 11 12 13 14 14 15 16 17 18 19 19 20 21 22
26 8 9 9 10 11 11 12 13 14 15 16 17 17 18 19 20 21 22 23 24 26
27 10 10 11 12 12 13 14 15 16 17 18 19 20 21 22 24 25 26 27 29 30
28 11 12 13 14 15 16 17 18 19 20 21 22 24 25 26 28 29 30 32 34 35
29 14 15 16 17 18 19 20 21 23 24 25 26 28 29 31 33 34 36 38 40 42
30 17 18 19 21 22 23 24 26 27 28 30 31 33 35 37 39 41 43 45 48 51
31 21 22 24 25 27 28 29 31 32 34 36 37 39 41 43 46 48 51 54 57 61
32 25 27 29 30 32 34 36 37 39 41 42 44 46 49 51 54 57 61 65 69 73
33 30 32 35 37 39 41 43 45 46 48 50 52 55 57 61 64 68 72 76 81 86
34 36 39 41 44 46 49 51 53 55 57 59 61 64 67 71 75 79 84 89 95 101
35 42 45 48 51 54 57 59 61 63 66 68 71 74 78 82 86 91 97 103 109 116
36 48 52 55 59 62 65 67 70 72 75 78 81 85 89 93 98 104 110 116 123 130
37 54 58 62 65 69 72 75 78 81 84 87 91 95 99 104 110 116 122 129 136 143
38 59 63 67 71 75 78 82 86 89 93 97 101 105 110 115 121 127 133 140 147 155
39 63 67 72 76 80 84 88 92 97 101 105 110 115 120 126 131 138 144 151 157 165
40 68 72 76 81 85 90 95 99 104 109 114 120 125 131 136 142 148 154 161 167 174
41 73 78 82 87 92 97 102 107 113 119 125 130 136 142 148 154 160 166 173 179 186
42 80 85 89 94 100 105 111 117 123 130 136 143 149 155 161 168 174 180 186 192 199
43 89 94 99 104 110 116 122 129 136 143 150 157 163 170 177 183 189 195 202 208 215
44 100 105 110 116 122 129 136 143 150 158 165 173 180 187 194 200 207 213 220 226 233
45 113 119 125 132 138 145 153 160 168 176 183 191 198 206 213 220 227 234 241 248 255
46 131 137 144 151 158 166 173 181 188 196 204 212 219 227 234 242 249 257 265 273 281
47 153 161 168 176 183 191 198 206 213 220 228 235 243 251 259 267 275 283 292 301 310
48 181 190 199 207 215 222 229 235 242 248 255 262 269 277 285 294 303 313 323 333 344
49 217 227 237 245 253 259 265 270 275 279 285 291 298 306 315 324 335 346 357 369 381
50 259 272 283 291 298 303 306 309 311 313 317 322 329 336 346 356 368 380 394 407 420
51 307 321 333 341 346 349 349 349 349 349 351 354 360 368 378 389 402 416 431 446 461
52 355 371 383 390 394 394 392 389 386 384 385 387 393 401 411 423 437 452 468 484 500
53 399 416 429 435 437 436 432 427 422 419 417 420 425 433 443 456 470 486 503 520 536
54 433 452 464 470 472 470 465 459 454 450 448 450 456 464 474 487 502 518 535 553 569
55 456 473 486 492 494 493 489 485 480 477 476 479 484 493 504 517 532 548 564 581 597
56 465 482 494 502 506 507 505 503 501 500 501 504 511 520 531 544 558 574 590 605 620
57 464 479 492 501 507 511 514 515 516 518 522 527 535 544 556 568 582 596 611 625 638
58 457 471 484 495 504 511 517 522 528 533 539 547 556 566 578 590 603 616 629 641 653
59 449 462 475 487 498 509 518 528 537 546 555 565 576 587 598 610 621 633 644 655 665
60 444 456 470 482 495 508 521 533 546 558 571 583 594 606 617 628 639 649 659 668 677
61 445 457 470 484 498 512 527 542 558 572 587 600 613 625 636 647 657 666 674 682 689
62 452 464 477 491 506 522 538 555 572 588 604 619 632 645 656 665 674 682 690 696 703
63 464 476 489 503 519 535 553 571 588 606 623 638 652 664 674 684 692 699 706 712 718
64 482 493 507 521 537 554 571 589 608 625 642 657 671 683 693 702 710 716 722 728 733
65 505 516 529 543 559 576 593 611 629 646 663 677 691 702 712 720 727 734 739 744 749
66 532 544 557 570 586 601 618 635 652 668 684 698 710 721 730 738 745 751 756 761 766
67 564 576 588 601 616 630 646 661 677 691 705 718 729 739 748 755 762 768 773 778 783
68 600 611 623 635 649 662 676 689 702 715 727 738 748 757 765 772 778 784 790 795 800
69 639 649 660 671 683 694 706 717 728 738 748 757 766 774 781 788 794 800 806 811 816
70 677 687 697 707 717 726 735 744 752 760 768 776 783 790 797 803 810 816 821 827 832
71 712 722 731 739 748 755 762 768 775 781 787 793 799 805 812 818 824 830 836 841 846
72 744 753 761 768 775 781 786 790 795 799 804 809 814 820 826 831 837 843 849 854 859
73 770 779 786 792 798 803 807 810 813 816 820 824 829 833 839 844 850 855 861 866 871
74 792 800 806 812 817 821 824 826 828 831 834 838 842 846 851 856 861 866 872 877 881
75 808 816 821 827 831 835 837 840 842 844 847 850 854 858 862 867 872 877 882 886 891
76 821 828 833 838 842 846 849 851 853 855 858 861 865 869 873 877 882 886 891 895 899
77 829 836 841 846 851 855 857 860 863 865 868 871 875 879 883 887 891 895 899 903 906
78 835 842 847 852 857 861 864 868 871 874 877 881 884 888 892 895 899 903 906 910 913
79 840 846 851 856 862 867 870 874 878 882 886 889 893 896 900 903 907 910 913 916 919
80 844 850 855 861 866 872 876 881 885 889 893 897 901 904 908 911 914 917 919 922 924
81 849 854 860 865 871 877 882 887 892 896 901 905 908 912 915 918 921 923 925 928 930
82 855 860 866 871 877 883 888 893 898 903 908 912 916 919 922 924 927 929 931 933 935
83 862 866 872 877 883 889 894 900 905 910 915 919 922 925 928 931 933 935 937 938 940
84 869 873 879 884 890 896 901 906 912 917 921 925 929 932 934 936 938 940 942 943 945
85 878 881 887 891 897 903 908 913 918 923 927 931 935 937 940 942 944 945 947 948 949
86 886 890 895 899 905 910 915 920 925 930 934 937 940 943 945 947 949 950 951 953 954
87 896 899 904 908 913 918 923 927 932 936 940 943 946 948 950 952 954 955 956 957 958
88 906 908 913 917 922 926 930 934 939 942 946 949 951 953 955 957 958 960 961 962 963
89 916 918 923 926 931 934 938 942 945 949 952 954 957 959 960 962 963 964 965 966 967
90 927 929 933 935 940 942 946 949 952 955 958 960 962 964 965 967 968 969 970 971 972

basically, here how the table works, if i measured 317 on ADC at 25deg centigrade, it means, RH is 50%RH (pa trace nalang po sa table). Thanks to sir zero for helping me with the coding.

I wanted to get the function of the curve. baka may nakaexperience na sainyo ng pag-gawa nito. i made use of datafit9 by oakdale engineering to generate those values. i use rational function interpolation.  kung pupwede using matlab kasi may curve fitting tool din naman yun.


NOTE: working na po yung actual humidity/temperature sensor node. and comparing the result to the electronic one, ok na.
one of this day, i will post a simple project using CM-R humidity sensor.

[insomartin]

that's a lot of data to look up... like 5kbytes...
i'm having problems with non/semi log, acceleration curves.
but that's me. i wanted an interpolated data to show like 82.01 % than 82&

[insomartin]

^subscribing

i also wanted to know how to make the formulas. thanks

[marcelino]

sorry, di pala "rational function interpolation" yung ginamit ko. cubic spline pala.

last night, natry ko sa matlab yung cubic spline. very simple.

this is the content of the data.txt. the resistance of the sensor at 25deg.

20   5000
30   1300
40   310
50   87
60   31
70   13
80   5.7
90   2.7

here's what I got, based on the HELP.



however, i wanted to get the values for each point [20, 21,22,23....]

Code:

load data.txt;
x=data(:,1);
y=data(:,2);
cs = csapi(x,y);
fnplt(cs);
hold on;
plot(x,y,'o');
fnplt(cs);
hold off;

alam nyo ba yung command ng pagkuha nung values for each point that i wanted? [21,22,23,24...]


inso, try mo yung sa oakdale engineering. trial version could do that. but then again, i want to get the equation...
heheh

[insomartin]

OO NGA!!! cubic spine thanks thanks! bwahahahaha tagal ko na nag cocoreldraw! yun nga tmang term dun!

pero paano makuha yung interpolated splines?  

[e.novacek]

alam nyo ba yung command ng pagkuha nung values for each point that i wanted? [21,22,23,24...]

Ans:

we can re-define the m-file as follows:

>> x = [ 20: 1 : 30 ]   % step count (+1 increment) starting from 20 ending to 30
                               % for abcissa
....
......

>> plot(x,y,'g+')         % plots the loaded *.txt interpolated values..

[zer0w1ng]

With 8 results of known values, a 8 equation/ 8 unknown equation could be deduced.
A 7th order rh(resistance) function could be expressed this way :
rh(x) = a*x^7 + b*x^6 + c*x^5 + d*x^4 + e*x^3 + f*x^2 + g*x + h

Solving using a math solver program (mathcad could be used with this):
rh(10000) = 20
rh(2330) = 30
rh(550) = 40
rh(155) = 50
rh(51) = 60
rh(20) = 70
rh(8.3) = 80
rh(4.1) = 90

21 = rh(r21)
22 = rh(r22)
23 = rh(r23)

Results to:
  a = -8.7832161783323077E-16
  b = +1.1524127623330964E-11
  c = -2.9156944955937892E-08
  d = +1.7788599438337333E-05
  e = -0.0031212946621026888
  f  = +0.18393856255949986
  g = -4.3067170086982713
  h        = +104.7756623408

So the rh values of 21,22, 23 and 24 are:
  r21      = +67.7210967157963
  r22      = +67.4078288641299
  r23      = +67.0915638552075
  r24      = +66.772188350739

[marcelino]

alam nyo ba yung command ng pagkuha nung values for each point that i wanted? [21,22,23,24...]

Ans:

we can re-define the m-file as follows:

>> x = [ 20: 1 : 30 ]   % step count (+1 increment) starting from 20 ending to 30
                               % for abcissa
....
......

>> plot(x,y,'g+')         % plots the loaded *.txt interpolated values..

thanks... pero di kaya ang result maging straight line from 20 to 30...?

With 8 results of known values, a 8 equation/ 8 unknown equation could be deduced.
A 7th order rh(resistance) function could be expressed this way :
rh(x) = a*x^7 + b*x^6 + c*x^5 + d*x^4 + e*x^3 + f*x^2 + g*x + h

Solving using a math solver program (mathcad could be used with this):
rh(10000) = 20
rh(2330) = 30
rh(550) = 40
rh(155) = 50
rh(51) = 60
rh(20) = 70
rh(8.3) = 80
rh(4.1) = 90

21 = rh(r21)
22 = rh(r22)
23 = rh(r23)

Results to:
  a = -8.7832161783323077E-16
  b = +1.1524127623330964E-11
  c = -2.9156944955937892E-08
  d = +1.7788599438337333E-05
  e = -0.0031212946621026888
  f  = +0.18393856255949986
  g = -4.3067170086982713
  h        = +104.7756623408

So the rh values of 21,22, 23 and 24 are:
  r21      = +67.7210967157963
  r22      = +67.4078288641299
  r23      = +67.0915638552075
  r24      = +66.772188350739

sige po, try ko to mamaya... download muna ng mathcad.

pero observation ko po... bakit magkakalapit yung resulted values. 67.72 - 67.40 - 67.01?

malaki po yung drop form 20 to 24(near 25), according sa graph.

[zer0w1ng]

As observed, RH of 70-80 has large drop fro 20ohms to 8.3ohms which is in very non-linear region.
So a resistance of 21 will be very near to 70%RH as per result of the curve matching.

[marcelino]

I just got how its going to be.

Ang issue ko naman ngayon ay kung alin ang gagamitin kong curve fitting.

May curve fitting tool pala ang MATLAB. nakalimutan ko...
after loading data as X and Y. i used the data on 15deg centigrade.

I tried to use a 7th degree polynomial fitting. which as of now gave me a very good impression on the plotted curve.
eto po.


eto naman po yung result ng equation at ibang details... di ko pa maintindihan.

Warnings during fitting:
Equation is badly conditioned. Remove repeated data points
         or try centering and scaling.
 
Linear model Poly7:
       f(x) = p1*x^7 + p2*x^6 + p3*x^5 + p4*x^4 + p5*x^3 +
                    p6*x^2 + p7*x + p8
Coefficients:
       p1 = -2.448e-008
       p2 =   1.11e-005
       p3 =   -0.002146
       p4 =      0.2297
       p5 =      -14.74
       p6 =       569.3
       p7 = -1.231e+004
       p8 =  1.159e+005

Goodness of fit:
  SSE: 2.953e-018
  R-square: 1
  Adjusted R-square: NaN
  RMSE: NaN

how does one uses the "goodness of fit?"

[e.novacek]

the result will be straight-line not unless the "ordinate" is linear with respect to x..

to illustrate what i mean, here's a simple program: (under matlab command prompt..)
...
>>
>> r = 0: pi/50: 5*pi ;
>>
>> s = 10*sin(2*pi*r);
>>
>> plot(r,s,’r*’);
>>
... ....

in the above example.. the range is extending from 0 to 5*pi however, there is a "post increment" between them where it numerically define as pi/50.. literally "r" is piece-wise linear..

non-linearity occurs in this case when "s" intervenes with respect to r..  

my suggestion for clarity of the plot:

 %=========================
 % include in m-file..
 %=========================
  ....
  .......
 xlabel(‘ state the name of x ’ );
 ylabel(‘ state the name of y ’ );
 title( ‘ title of the program ‘ );
 grid on ;

i suggest to use the mesh command to visualize clearly the interpolation.. in accordance with the graph.. the counting is by 10's and that is why na-omit ung 21,22,23..

so kung gusto natin i-include ung 21, 22, 23.. then we should changed the incrementing of numbers.. let's say:

>> T = 20: 1: 30 ..
>> ..

then magbabago na ung counting ng plot from 20 ~ 90 wherein ang increment by 1's na.. so dapat baguhin din ung counting (or interpolation) ng ordinate-axes.. (which is dependent on the equation usually y-axes) then the effects of counting will be plotted in graph..
 
i think much better kung excel ang gamitin sa plotting when "floating points" are involved..

[marcelino]

the result will be straight-line not unless the "ordinate" is linear with respect to x..

to illustrate what i mean, here's a simple program: (under matlab command prompt..)
...
>>
>> r = 0: pi/50: 5*pi ;
>>
>> s = 10*sin(2*pi*r);
>>
>> plot(r,s,’r*’);
>>
... ....

in the above example.. the range is extending from 0 to 5*pi however, there is a "post increment" between them where it numerically define as pi/50.. literally "r" is piece-wise linear..

non-linearity occurs in this case when "s" intervenes with respect to r..  

my suggestion for clarity of the plot:

 %=========================
 % include in m-file..
 %=========================
  ....
  .......
 xlabel(‘ state the name of x ’ );
 ylabel(‘ state the name of y ’ );
 title( ‘ title of the program ‘ );
 grid on ;

i suggest to use the mesh command to visualize clearly the interpolation.. in accordance with the graph.. the counting is by 10's and that is why na-omit ung 21,22,23..

so kung gusto natin i-include ung 21, 22, 23.. then we should changed the incrementing of numbers.. let's say:

>> T = 20: 1: 30 ..
>> ..

then magbabago na ung counting ng plot from 20 ~ 90 wherein ang increment by 1's na.. so dapat baguhin din ung counting (or interpolation) ng ordinate-axes.. (which is dependent on the equation usually y-axes) then the effects of counting will be plotted in graph..
 
i think much better kung excel ang gamitin sa plotting when "floating points" are involved..

di ko pa din masundan sis...

eto yung data ko:

Code:

20 10000
30 2330
40 550
50 155
60 51
70 20
80 8.3
90 4.1

So after loading that data:

Code:

x =

    20
    30
    40
    50
    60
    70
    80
    90

and

Code:

y =

  1.0e+004 *

    1.0000
    0.2330
    0.0550
    0.0155
    0.0051
    0.0020
    0.0008
    0.0004

so sunod is
>>plot(x,y,'r*');

papano ko ilalagay dito yung increment of 1. nagawa ko nga using x = [20:1:90], pero di na magkapareho ang size ng dalawang array...

di ko talga masundan. heheh

[marcelino]

hehehe

i was able to use cubic spline in the curve fitting tool of matlab (cftool), pero walang naiprovide na equation. at yung result ay pareho lang sa paggawa ng spline function.
eto nalang ang ginawa ko:

%Generated Interpolation using Cubic Spline Using Matlab
%yy=spline(x,y,xx)
%example: data at 25degree
x=[20,30,40,50,60,70,80,90]
y=[5000,1300,310,87,31,13,5.7,2.7]
yy25deg = spline(x,y,[20:1:30]);

the result should be (transposed, yy25deg' ):

Code:

5000
4443.32965
3934.197305
3470.332076
3049.463074
2669.319408
2327.630189
2022.124529
1750.531537
1510.580324
1300
1116.519676
957.8684632
821.7754711
705.9698105
608.1805921
526.1369263
457.5679237
400.2026947
351.77035
310
272.9486447
239.9848421
210.8050395
185.1056842
162.5832237
142.9341053
125.8547763
111.0416842
98.19127632
87
77.19574474
68.63216842
61.19437105
54.76745263
49.23651316
44.48665263
40.40297105
36.87056842
33.77454474
31
28.44937632
26.09448421
23.92447632
21.92850526
20.09572368
18.41528421
16.87633947
15.46804211
14.17954474
13
11.91945
10.93149474
10.03062368
9.211326316
8.468092105
7.795410526
7.187771053
6.639663158
6.145576316
5.7
5.297423684
4.932336842
4.599228947
4.292589474
4.006907895
3.736673684
3.476376316
3.220505263
2.96355
2.7

yan na yata ang sinasabi ni sis e.

pero ganun pa din walang equation. pero yan na nga yata yan. kaya nga gumagamit ng MATLAB para mapadali ang calculations.

[zer0w1ng]

pero ganun pa din walang equation.

The solution I gave from my previous post resulted already in an equation.
It is a 7th order polynomial:

    rh(x) =  a*x^7 + b*x^6 + c*x^5 + d*x^4 + e*x^3 + f*x^2 + g*x + h

and the coefficients a to h is as follows:

  a = -8.7832161783323077E-16
  b = +1.1524127623330964E-11
  c = -2.9156944955937892E-08
  d = +1.7788599438337333E-05
  e = -0.0031212946621026888
  f  = +0.18393856255949986
  g = -4.3067170086982713
  h = +104.7756623408

I guess I have not explained well my post before.

[piona]

Wala po bang ideal Curve Equation na ibinigay?

Kasi pagginawa na ito kailangan e calibrate yong device response, againts sa ideal curve.
 

[marcelino]

The solution I gave from my previous post resulted already in an equation.
It is a 7th order polynomial:

    rh(x) =  a*x^7 + b*x^6 + c*x^5 + d*x^4 + e*x^3 + f*x^2 + g*x + h

and the coefficients a to h is as follows:

  a = -8.7832161783323077E-16
  b = +1.1524127623330964E-11
  c = -2.9156944955937892E-08
  d = +1.7788599438337333E-05
  e = -0.0031212946621026888
  f  = +0.18393856255949986
  g = -4.3067170086982713
  h = +104.7756623408

I guess I have not explained well my post before.

sir, yung ginawa ko po ay interpolation lamang. its called "cubic spline".

may issue po ako sa 7th order polynomial curve fitting. di naman po nagfit.

example:
a = -8.7832161783323077E-16
  b = +1.1524127623330964E-11
  c = -2.9156944955937892E-08
  d = +1.7788599438337333E-05
  e = -0.0031212946621026888
  f  = +0.18393856255949986
  g = -4.3067170086982713
  h = +104.7756623408

x = 50

f(x) = [should be 115@15deg] but it resulted to  61.3710.

ginawa ko na din po yung 7th order polynomial using MATLAB.
eto po ang result for the 15degree centigrade.

Code:

f(x) = p1*x^7 + p2*x^6 + p3*x^5 + p4*x^4 + p5*x^3 +
                    p6*x^2 + p7*x + p8
       where x is normalized by mean 55 and std 24.49
Coefficients:
       p1 =      -129.5
       p2 =       361.9
       p3 =      -334.1
       p4 =       252.5
       p5 =      -374.3
       p6 =       383.6
       p7 =      -238.6
       p8 =       86.55

di din tumatama. whereas dun po sa cubic spline interpolation, sakto po.

Wala po bang ideal Curve Equation na ibinigay?

Kasi pagginawa na ito kailangan e calibrate yong device response, againts sa ideal curve.
 

wala eh. eto yung datasheet: http://i01.yizimg.com/upload/3381/20061221103057890280486.pdf

nagvavary kasi ang characteristic ng humidity sensor depende sa temperature.

[zer0w1ng]

x = 50
f(x) = [should be 115@15deg] but it resulted to  61.3710.
nagvavary kasi ang characteristic ng humidity sensor depende sa temperature.

Baliktad ang gusto mo, the function is with resistance as the independent variable and the dependent is RH.

    RH = f(resistance)

So r(50) = values between 60-70%RH, where 50 is the resistance.
Thus the equation correctly gives 61.3710%RH as the humidity.

This is way, having the independent variable as the resistance for this is what will be read by the MCU.
Passing to the function will give you the right RH.
If you want the reverse, we can re-compute again for the function.

[z22]

Halos linear yung mga curve sa graph ni Martin Mystère kung log yung scale ng resistance.
Paano kung kunin mo na lang yung log(resistance), tapos 2D linear interpolation, mas simple yun.

[insomartin]

it's log sir... 

[insomartin]

OT:
lagi na lang OT post ko