用C++做資料分析 | CERN ROOT 教學[03] - 函數與函數擬合


簡介

在上一篇文章中我們介紹了基本資料分析的操作方法。在這篇文章中我們會詳細介紹在資料分析中常見的一個功能:函數擬合

函數

ROOT的TF1可以用來表示任意的一元函數。他的自變數的名稱永遠為x。而且你可以在裡面使用任意的C/C++函數(包括ROOT自帶的TMath函數庫)。

例如創造一個範圍在[-10, 10]之間的sinc函數。

root [0] f1 = new TF1("f1", "sin(x)/x", -10, 10);
root [1] f1->Draw();
sin(x)/x

因為你可以用任何 C/C++的函數。所以也可以自己定義一個std::function或是lambda function後使用她。值得注意的是因為C++17支援Class template type deduction而且ROOT自動匯入std的命名空間。所以原本在C++11宣告一個std::function的複雜語句auto my_func = std::function<double(double)>([](double){...})可以簡化成my_func = function([](double){...})

root [2] my_sinc = function([](double x){ return sin(x)/x; });
root [3] f1 = new TF1("f1", "my_sinc(x)", -10, 10);
root [4] f1->Draw();
sin(x)/x

說到這裡熟悉C++的讀者應該會想到:那具有function call operator的類別是否也可以當作函數使用呢?答案是可以的。我在這裡就不多做說明,不過應該很容易成功。

有了函數,還需要能夠設定函數的各種參數。像是想要表達一個任意的線性函數y = ax+b的時候需要能夠指定固定的a與b。我們可以在函數的表達式中加入[0]來代表第0個參數、[1]代表第一個,以此類推。並在繪製圖形前用TF1::SetParameters來設定數值(其型別固定為double)

root [5] f1 = new TF1("f1", "[0]*x+[1]", -10, 10);
root [6] f1->SetParameters(-1, -5);
root [7] f1->Draw();
y = ax+b,a = -1, b = -5

你也可以使用變數名稱來指定變數。這種方法比較好寫也比較容易除錯。但不像使用索引一樣使用SetParameters一次指定一堆變數。

root [8] f1 = new TF1("f1", "exp(x+[a])", -10, 10);
root [9] f1->SetParameter("a", 2);
root [10] f1->Draw();
exp(x+a)

當需要取得變數的數值的時候可以使用GetParameter來取得函數。

root [11] f1->GetParameter("a")
(double) 2.0000000

函數擬合

講完函數之後終於可以討論如何使用函數擬合了!首先。我們先生出一個常態分佈。

root [0] h1 = new TH1F("h1", "histogram", 60, -3, 3);
root [1] h1->FillRandom("gaus", 100000);
root [2] h1->Draw();
常態分佈

並且如前一篇文章所述,使用h1->Fit("gaus");將高斯分佈函數擬合上去。

root [4] h1->Fit("gaus");
 FCN=58.2284 FROM MIGRAD    STATUS=CONVERGED      62 CALLS          63 TOTAL
                     EDM=8.58021e-13    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  Constant     4.00036e+03   1.57642e+01   4.75510e-02  -9.85101e-09
   2  Mean        -1.17823e-03   3.20398e-03   1.20393e-05  -2.62278e-05
   3  Sigma        9.99354e-01   2.41238e-03   2.44449e-06  -1.71898e-03
將高斯分佈擬合至高斯分佈上

之後我們可以使用GetFunction取得擬合後函數本身並取得函數的所有參數。在此函數的名稱為你所使用的分佈,而參數的名稱在上面擬合後的輸出結果中。

root [5] f1 = h1->GetFunction("gaus");
root [6] f1->GetParameter("Constant")
(double) 4000.3641
root [7] f1->GetParameter("Sigma")
(double) 0.99935417

擬合自定義函數

當然我們不只能擬合到預先定義好的函數上,也可以擬合到自定義的函數。不過有兩個前提:要擬合的函數必須帶有參數,且函數必須是determinstic的。所以ROOT可以幫你擬合[scale]*x+[bias]+exp(x*[exp_scale])但沒辦法擬合(rand()%1000)*x*[scale]

要擬合自定義函數需要先定義一個TF1函數。然後呼叫TH1::Fit將函數擬合上去。

root [7] f2 = new TF1("f2", "[scale]*sin(x)/x+[bias]");
root [8] h1->Fit(f2);
 FCN=14310 FROM MIGRAD    STATUS=CONVERGED      30 CALLS          31 TOTAL
                     EDM=7.17638e-19    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  bias        -3.32813e+02   3.54589e+00   1.38923e-01   2.29136e-10
   2  scale        2.85750e+03   1.11800e+01   4.38015e-01  -4.15280e-11
擬合自定義函數

小結

我們可以在ROOT中輕松的建立數學函數,並可以將任意函數擬合到分佈上。函數擬合是在資料分析時相當常用的功能。

#CERN, ROOT







你可能感興趣的文章

菜比八寫後端(1) - MySQL簡介

菜比八寫後端(1) - MySQL簡介

關於 JavaScript 的提升 (Hoisting)

關於 JavaScript 的提升 (Hoisting)

React-[表單篇]-表單開發|手工版本

React-[表單篇]-表單開發|手工版本






留言討論