最後の場合を、TH1を使って簡単に求める方法を考える。
鍵はTH1::IntegralAndError().
以下のコードに例を示す。 適当なweight wを与えた10イベントに対して、method 1は√( Σ(w_i)^2 )の顕な計算によって誤差を求めたもの、 method 2はTH1を使った方法。 共に同じ値を返すので、method 2で正しく誤差が求められていることが確かめられる。
void IntegralError() { const Int_t N = 10; const Double_t w[N] = { 0.9, 0.9, 0.9, 0.9, 0.9, 1.0, 1.0, 1.0, 1.0, 1.0 }; {/// method 1 Double_t Integral = 0; Double_t eIntegral = 0; for(Int_t i = 0 ; i < N ; ++i){ Integral += w[i]; eIntegral += pow( w[i], 2 ); } eIntegral = sqrt(eIntegral); ///// result cout << "method 1" << endl; cout << "Integral = " << Integral <<" +/- "<< eIntegral << endl; } {/// method 2 TH1D *h = new TH1D("h","h",1,0,1); h->Sumw2(); for(Int_t i = 0 ; i < N ; ++i) h->Fill( 0., w[i] );///// 0のbinにweightをかけて詰める Double_t eIntegral; Double_t Integral = h->IntegralAndError(0,1+1,eIntegral);/////積分はunderflow(bin=0),overflow(bin=2)も含めるようにしている(一応) ///// result cout << "method 2" << endl; cout << "Integral = " << Integral <<" +/- "<< eIntegral << endl; } }
例えばtreeからある条件を満たすweight付きのイベント数およびその統計誤差を求めたければ、
const TCut cut = "...";///// カット条件 const TCut weight = "...";///// weight TH1D *h = new TH1D("h","",1,0,1); h->Sumw2(); tree->Project("h","val",cut*weight);///// val は適当な変数名(treeが持っているbranch名の1つ) Double_t eIntegral; Double_t Integral = h->IntegralAndError(0,2,eIntegral); cout << " Integral = " << Integral << " +/- " << eIntegral << " ( " << eIntegral/Integral *100 << " % )" << endl;
のようにすれば良い。