ROOT 分析作业

查看练习示例脚本

    

     

      
     
     

     

     

      

       

        

         

          

          

          

          

           

            

            

            
>cdroot_exercise/
           
           
>ls
           
-rw-r--r--1user group1792Aug1619:45 exercise0.C
           
-rw-r--r--1user group969Aug1619:49 exercise1.C
           
-rw-r--r--1user group2848Aug1619:40 exercise2.C
          
         
        
       
      
      

      

     
    

1，生成数据样本

1.1 定义分布函数

Lorentzian Peak function:

    

     

      
     
     

     

     

      

       

        

         

          

           x
          
          

          

          

           

            

            

            
// Quadratic background function
           
           
Double_tbackground(Double_t*x,Double_t*par) {
           
returnpar[0]+par[1]*x[0]+par[2]*x[0]*x[0];
           
}
           

           
// Lorentzian Peak function
           
Double_tlorentzianPeak(Double_t*x,Double_t*par) {
           
return(0.5*par[0]*par[1]/TMath::Pi())/TMath::Max(1.e-10,
           
(x[0]-par[2])*(x[0]-par[2])+.25*par[1]*par[1]);
           
}
           

           
// Sum of background and peak function
           
Double_tFunction(Double_t*x,Double_t*par) {
           
returnbackground(x,par)+lorentzianPeak(x,&par[3]);
           
}
          
         
        
       
      
      

      

     
    

1.2 生成数据样本，填tree，写入root文件

    

     

      
     
     

     

     

      

       

        

         

          

           
x
          
          

          

          

           
// create the file, the Tree and a few branches
           

            
TFile*f=TFile::Open("tree.root","recreate");
           
           
TTree*tree=newTTree("tree","treelibrated tree");
           
Float_tx,y;
           

           
// create a branch with energy
           
tree->Branch("X",&x);
           
tree->Branch("Y",&y);
           

           
// create a TF1 with the range from 0 to 3 and 6 parameters
           
TF1*Fcnx=newTF1("Fcnx",Function,0,3,6);
           
TF1*Fcny=newTF1("Fcny",Function,0,3,6);
           

           
// Fix the Parameters of function
           
Fcnx->SetParameters(-1,45,-13.3,13.8,0.2,1);
           
Fcny->SetParameters(-1,45,-13.3,13.8,0.1,1.5);
           

           
// fill some events with random numbers
           
Int_tnevent=1000;
           
for(Int_tiev=0;iev<nevent;iev++) {
           
x=Fcnx->GetRandom();
           
y=Fcny->GetRandom();
           
tree->Fill();// fill the tree with the current event
           
}
           

           
// save the Tree header. The file will be automatically closed
           
// when going out of the function scope
           

            

            

            
tree.Write();
           
          
         
        
       
      
      

      

     
    

操作练习

    

     

      
     
     

     

     

      

       

        

         

          

           
xxxxxxxxxx
          
          

          

          

           

            

            

            
> root-lexercise0.C
           
           
root [0]
           
Processing exercise0.C...
           
root [1] .q
           
# 查看生成的 root 文件
           
>ls-l
           
-rw-r--r--1user group1792Aug1619:45 exercise0.C
           
-rw-r--r--1user group969Aug1619:49 exercise1.C
           
-rw-r--r--1user group2848Aug1619:40 exercise2.C
           
-rw-r--r--1user group76796Aug1620:03 tree.root
           

           
# 打开 root 文件, 通过 TBrowser 查看内容
           
> root-ltree.root
           
root [0] TBrowser a
           
# 在弹出窗口左侧文件列表中 "ROOT Files" 目录下找到 tree.root 文件，
           
# 双击 root 文件查看tree "tree"，双击 tree 查看 Branch "X" 和 “Y”，
           
# 双击 Branch “X” 可以在右侧窗口中看到 “X” 的直方图。
          
         
        
       
      
      

      

     
    

2. 读取root文件，填直方图，画出直方图

    

     

      
     
     

     

     

      

       

        

         

          

           
x
          
          

          

          

           

            
TFile*f=TFile::Open("tree.root");
           
           
TTree*t1=(TTree*)f->Get("tree");
           

           
Float_tx,y;
           
t1->SetBranchAddress("px",&px);
           
t1->SetBranchAddress("py",&py);
           

           
// create two histograms
           
TH1F*hx=newTH1F("hx","x distribution",100,-3,3);
           
TH2F*hxy=newTH2F("hxy","y vs x",30,-3,3,30,-3,3);
           

           
// read all entries and fill the histograms
           
Long64_tnentries=t1->GetEntries();
           
for(Long64_ti=0;i<nentries;i++) {
           
t1->GetEntry(i);
           
hx->Fill(x);
           
hxy->Fill(x,y);
           
}
           

           
// draw the histograms
           
TCanvas*c1=newTCanvas();
           
hx->Draw("E0");
           
c1->Print("hx.png")// eps, ps, jpg
           

           
TCanvas*c2=newTCanvas();
           
hxy->Draw("colz");// BOX, SCAT, ARR
           

            
c1->Print("hxy.png")
           
          
         
        
       
      
      

      

     
    

操作练习

    

     

      
     
     

     

     

      

       

        

         

          

           
xxxxxxxxxx
          
          

          

          

           

            

            

            
> root-lexercise1.C
           
           
root [0]
           
# 可以看到生成了一个一维直方图以及二维直方图
           
# 尝试改变画图选项，得到不同风格的直方图
          
         
        
       
      
      

      

     
    

单击菜单栏中 View 下拉菜单选中 Editor 选项可以调出编辑面版，可以手动编辑图片属性，如坐标轴标题，添加网格，直方图颜色等。尝试用鼠标单击图像不同位置，如直方图，坐标轴，标题处，观察编辑面板编辑选项的变化。尝试手动编辑图像属性。

3. 拟合直方图

    

     

      
     
     

     

     

      

       

        

         

          

           
x
          
          

          

          

           
// Read the Tree
           
TFile*f=TFile::Open("tree.root");
           
TTree*t1=(TTree*)f->Get("tree");
           

           
TCanvas*c1=newTCanvas("c1","Fitting Demo",10,10,700,500);
           
c1->SetFillColor(33);
           
c1->SetFrameFillColor(41);
           
c1->SetGrid();
           

           
// Create one histogram
           
TH1F*histo=newTH1F("histo","Lorentzian Peak on Quadratic Background",60,0,3);
           
histo->SetMarkerStyle(21);
           
histo->SetMarkerSize(0.8);
           
histo->SetStats(0);
           

           
// Fill the histograms
           
t1->Project("histo","X");
           

           

            
// create a TF1 with the range from 0 to 3 and 6 parameters
           
           
TF1*fitFcn=newTF1("fitFcn",fitFunction,0,3,6);
           
fitFcn->SetNpx(500);
           
fitFcn->SetLineWidth(4);
           
fitFcn->SetLineColor(kMagenta);
           

           
// first try without starting values for the parameters
           
// This defaults to 1 for each param.
           
// this results in an ok fit for the polynomial function
           
// however the non-linear part (lorenzian) does not
           
// respond well.
           
fitFcn->SetParameters(1,1,1,1,1,1);
           
histo->Fit("fitFcn","0");
           

           
// second try: set start values for some parameters
           
fitFcn->SetParameter(4,0.2);// width
           
fitFcn->SetParameter(5,1);// peak
           

           
histo->Fit("fitFcn","V+","ep");
           
// improve the picture:
           
TF1*backFcn=newTF1("backFcn",background,0,3,3);
           
backFcn->SetLineColor(kRed);
           
TF1*signalFcn=newTF1("signalFcn",lorentzianPeak,0,3,3);
           
signalFcn->SetLineColor(kBlue);
           
signalFcn->SetNpx(500);
           
Double_tpar[6];
           

           
// writes the fit results into the par array
           
fitFcn->GetParameters(par);
           

           
backFcn->SetParameters(par);
           
backFcn->Draw("same");
           

           
signalFcn->SetParameters(&par[3]);
           
signalFcn->Draw("same");
           

           
// draw the legend
           
TLegend*legend=newTLegend(0.6,0.65,0.88,0.85);
           
legend->SetTextFont(72);
           
legend->SetTextSize(0.04);
           
legend->AddEntry(histo,"Data","lpe");
           
legend->AddEntry(backFcn,"Background fit","l");
           
legend->AddEntry(signalFcn,"Signal fit","l");
           
legend->AddEntry(fitFcn,"Global Fit","l");
           

            

            

            
legend->Draw();
           
          
         
        
       
      
      

      

     
    

操作练习

    

     

      
     
     

     

     

      

       

        

         

          

           
xxxxxxxxxx
          
          

          

          

           

            

            

            
> root-lexercise1.C
           
          
         
        
       
      
      

      

     
    

得到拟合结果：

可以在输出 log 里面找到函数参数的拟合值：

    

     

      
     
     

     

     

      

       

        

         

          

           
xxxxxxxxxx
          
          

          

          

           

            

            

            
NO. NAME VALUE ERROR NEGATIVE POSITIVE
           
           
1p0-6.98862e+001.03331e+00
           
2p12.36951e+025.46818e+00
           
3p2-6.93593e+012.12246e+00
           
4p38.09132e+015.59809e+00
           
5p42.05657e-011.77811e-02
           
6p51.01006e+005.72858e-03
          
         
        
       
      
      

      

     
    

也可以通过以下方式拿到单个参数值及其误差：

    

     

      
     
     

     

     

      

       

        

         

          

           
xxxxxxxxxx
          
          

          

          

           

            

            

            
# Get Associated Function
           
           
root[0] TF1 *fit=histo->GetFunction("fitFcn");
           
# value of the first parameter
           
root[1] Double_t p1=fit->GetParameter(0);
           
# error of the first parameter
           
root[2] Double_t e1=fit->GetParError(0);
           
# 同时我们还可以计算拟合的 chisquare 值
           
root[3] Double_t chi2=fit->GetChisquare();
          
         
        
       
      
      

      

     
    

附加练习

尝试将本例中的朗道分布替换为高斯分布，完成 root 文件生成，root 文件读入，画直方图，拟合直方图。