Page MenuHomeHEPForge
Files F19301858

D55.1759174522.diff
No OneTemporary
Actions

Size
11 KB
Referenced Files
None
Subscribers
None

D55.1759174522.diff
View Options

diff --git a/inc/LauCalcChiSq.hh b/inc/LauCalcChiSq.hh
--- a/inc/LauCalcChiSq.hh
+++ b/inc/LauCalcChiSq.hh
@@ -23,8 +23,12 @@
*/
#include "TH2Poly.h"
+#include "TGraph.h"
#include "TString.h"
+#include "ROOT/RDataFrame.hxx"
+#include "ROOT/RResultPtr.hxx"
+
#include <vector>
/*! \file LauCalcChiSq.hh
@@ -115,6 +119,8 @@
TString treeName1_;
//! Name of the high stats data tree
TString treeName2_;
+ //! Experiment number in the toy file
+ TString iExpt_;
//! Name of the x-coordinate branch in tree 1
TString xName1_;
@@ -124,6 +130,10 @@
TString yName1_;
//! Name of the y-coordinate branch in tree 2
TString yName2_;
+ //! Title of the x-coordinate for drawing
+ TString xTitle_;
+ //! Title of the y-coordinate for drawing
+ TString yTitle_;
//! The minimum bin content
Float_t minContent_;
@@ -140,7 +150,13 @@
TH2Poly* chiSqHisto_;
//! Histogram (constructed from template) filled with signed chisq of tree1 vs tree2
TH2Poly* chiSqSignedHisto_;
-
+ //! Graph for data
+ ROOT::RDF::RResultPtr<TGraph> g_data_smartptr;
+ TGraph* g_data_;
+ //! Graph for toy
+ ROOT::RDF::RResultPtr<TGraph> g_toy_smartptr;
+ TGraph* g_toy_;
+
//! Minimum x coordinate of histograms
Float_t xMin_;
//! Maximum x coordinate of histograms
diff --git a/src/LauCalcChiSq.cc b/src/LauCalcChiSq.cc
--- a/src/LauCalcChiSq.cc
+++ b/src/LauCalcChiSq.cc
@@ -192,9 +192,14 @@
histo1_->SetTitleSize(0.06,"y");
histo1_->SetLabelSize(0.05,"x");
histo1_->SetLabelSize(0.05,"y");
- histo1_->SetXTitle(xName1_);
- histo1_->SetYTitle(yName1_);
- histo1_->Draw("colz");
+ histo1_->SetMarkerColor(0);
+ g_data_->SetMarkerStyle(20);
+ g_data_->SetMarkerSize(.5);
+ g_data_->Draw("ap");
+ g_data_->GetXaxis()->SetTitle(xTitle_);
+ g_data_->GetYaxis()->SetTitle(yTitle_);
+ //histo1_->Draw("colztextsame");
+ histo1_->Draw("colzsame");
can.SaveAs("data.pdf");
histo2_->SetLabelFont(62,"x");
@@ -207,10 +212,20 @@
histo2_->SetLabelSize(0.05,"y");
histo2_->SetXTitle(xName1_);
histo2_->SetYTitle(yName1_);
- histo2_->Draw("colz");
- can.SaveAs("toy.pdf");
-
- if(-1.*pullHisto_->GetMinimum() > pullHisto_->GetMaximum()) pullHisto_->SetMaximum(-1.*pullHisto_->GetMinimum());
+ histo2_->SetMarkerColor(0);
+ g_toy_->SetMarkerStyle(20);
+ g_toy_->SetMarkerSize(.5);
+ g_toy_->Draw("ap");
+ g_toy_->GetXaxis()->SetTitle(xTitle_);
+ g_toy_->GetYaxis()->SetTitle(yTitle_);
+ //histo2_->Draw("colztextsame");
+ histo2_->Draw("colzsame");
+ can.SaveAs("toy.png");
+
+ if(-1.*pullHisto_->GetMinimum() > pullHisto_->GetMaximum()) {
+ pullHisto_->SetMinimum(pullHisto_->GetMinimum()); // Required to workaround a bug in TH2Poly where, if SetMinimum is not called, negative entries are not drawn by "Draw("text")"
+ pullHisto_->SetMaximum(-1.*pullHisto_->GetMinimum());
+ }
else pullHisto_->SetMinimum(-1.*pullHisto_->GetMaximum());
pullHisto_->SetLabelFont(62,"x");
@@ -221,10 +236,20 @@
pullHisto_->SetTitleSize(0.06,"y");
pullHisto_->SetLabelSize(0.05,"x");
pullHisto_->SetLabelSize(0.05,"y");
- pullHisto_->SetXTitle(xName1_);
- pullHisto_->SetYTitle(yName1_);
- pullHisto_->Draw("colz");
+ pullHisto_->SetMarkerColor(0);
+ pullHisto_->SetMinimum(-9.);
+ pullHisto_->SetMaximum(9.);
+ g_data_->SetTitle("(data - toy) / error");
+ g_data_->Draw("ap");
+ gStyle->SetPaintTextFormat(".2f");
+ pullHisto_->SetLineWidth(1);
+ pullHisto_->SetLineColor(kBlack);
+ pullHisto_->SetLineStyle(1);
+ //pullHisto_->Draw("colztextsame");
+ pullHisto_->Draw("colzsame");
+ g_data_->Draw("psame");
can.SaveAs("pull.pdf");
+ can.SaveAs("pull.C");
chiSqHisto_->SetLabelFont(62,"x");
chiSqHisto_->SetLabelFont(62,"y");
@@ -234,9 +259,10 @@
chiSqHisto_->SetTitleSize(0.06,"y");
chiSqHisto_->SetLabelSize(0.05,"x");
chiSqHisto_->SetLabelSize(0.05,"y");
- chiSqHisto_->SetXTitle(xName1_);
- chiSqHisto_->SetYTitle(yName1_);
- chiSqHisto_->Draw("colz");
+ chiSqHisto_->SetMarkerColor(0);
+ g_data_->Draw("ap");
+ //chiSqHisto_->Draw("colztextsame");
+ chiSqHisto_->Draw("colzsame");
can.SaveAs("chiSq.pdf");
chiSqSignedHisto_->SetLabelFont(62,"x");
@@ -247,9 +273,12 @@
chiSqSignedHisto_->SetTitleSize(0.06,"y");
chiSqSignedHisto_->SetLabelSize(0.05,"x");
chiSqSignedHisto_->SetLabelSize(0.05,"y");
- chiSqSignedHisto_->SetXTitle(xName1_);
- chiSqSignedHisto_->SetYTitle(yName1_);
- chiSqSignedHisto_->Draw("colz");
+ chiSqSignedHisto_->SetXTitle(xTitle_);
+ chiSqSignedHisto_->SetYTitle(yTitle_);
+ chiSqSignedHisto_->SetMarkerColor(0);
+ g_data_->Draw("ap");
+ //chiSqSignedHisto_->Draw("colztextsame");
+ chiSqSignedHisto_->Draw("colzsame");
can.SaveAs("chiSqSigned.pdf");
}
@@ -257,32 +286,26 @@
void LauCalcChiSq::initialiseHistos()
{
- // Open the input control file:
- // Low_stat_file_name Low_stat_tree_name Low_stat_x_axis_name Low_stat_y_axis_name
- // High_stat_file_name High_stat_tree_name High_stat_x_axis_name High_stat_y_axis_name
- // Min_bin_content N_free_params Low/high_stat_histo_ratio xMin xMax yMin yMax
+ // Open the input control file:
+ // Low_stat_file_name Low_stat_tree_name iExpt Low_stat_x_axis_name Low_stat_y_axis_name
+ // High_stat_file_name High_stat_tree_name High_stat_x_axis_name High_stat_y_axis_name
+ // Min_bin_content N_free_params Low/high_stat_histo_ratio xMin xMax yMin yMax
std::ifstream getData(inputFileName_.Data());
// get the info on the low stat histo
- getData >> fileName1_ >> treeName1_ >> xName1_ >> yName1_;
+ getData >> fileName1_ >> treeName1_ >> iExpt_ >> xName1_ >> yName1_ >> xTitle_ >> yTitle_;
+
+ // Hack to introduce spaces in axis titles
+ xTitle_=xTitle_.ReplaceAll("SPACE"," ");
+ yTitle_=yTitle_.ReplaceAll("SPACE"," ");
+
if (!getData.good()) {
std::cerr<<"Error. Could not read first line of the input file "<<inputFileName_<<std::endl;
gSystem->Exit(EXIT_FAILURE);
}
- // open the file that contains the low stat histogram
- TFile * file1 = TFile::Open(fileName1_.Data(), "read");
-
- if (file1 == 0) {gSystem->Exit(EXIT_FAILURE);}
-
- // retrieve the low stat histogram
- TTree* tree1 = dynamic_cast<TTree*>(file1->Get(treeName1_.Data()));
-
- if (tree1 == 0) {
- std::cerr<<"Error. Could not find the tree "<<treeName1_
- <<" in the file "<<fileName1_<<std::endl;
- gSystem->Exit(EXIT_FAILURE);
- }
+ // Open the low stat TTree (we'll be interested in one experiment only)
+ ROOT::RDataFrame d1(treeName1_.Data(), fileName1_.Data());
// get the info on the high stat histogram
getData >> fileName2_ >> treeName2_ >> xName2_ >> yName2_;
@@ -293,24 +316,13 @@
// if both histograms are in the same file then just retrieve the
// high stat histogram from the file we already have open
- TFile * file2(0);
- TTree* tree2(0);
+ ROOT::RDataFrame d2=0;
if ( fileName2_ == fileName1_ ) {
- tree2 = dynamic_cast<TTree*>(file1->Get(treeName2_.Data()));
+ d2 = d1;
}
// otherwise open the other file and retrieve the high stat histogram from there
else {
- file2 = TFile::Open(fileName2_.Data(), "read");
-
- if (file2 == 0) {gSystem->Exit(EXIT_FAILURE);}
-
- tree2 = dynamic_cast<TTree*>(file2->Get(treeName2_.Data()));
- }
-
- if (tree2 == 0) {
- std::cerr<<"Error. Could not find the tree "<<treeName2_
- <<" in the file "<<fileName2_<<std::endl;
- gSystem->Exit(EXIT_FAILURE);
+ d2 = ROOT::RDataFrame(treeName2_.Data(), fileName2_.Data());
}
// get the info on the minimum content, number of parameters and scalefactor
@@ -330,72 +342,50 @@
// close the text file
getData.close();
- std::cout<<"Using the files and trees: "<<fileName1_<<"; "<<treeName1_
- <<" and "<<fileName2_<<"; "<<treeName2_<<std::endl;
- std::cout<<"Relative scaling factor histo1/histo2 set to "<<scaleFactor_<<std::endl;
- std::cout<<"Minimum bin content is "<<minContent_<<std::endl;
- std::cout<<"Number of free parameters is "<<nParams_<<std::endl;
-
- Double_t x;
- Double_t y;
-
- tree1->SetBranchAddress(xName1_.Data(),&x);
- tree1->SetBranchAddress(yName1_.Data(),&y);
-
- Int_t nEntries = tree1->GetEntries();
+ std::cout <<"Using the files and trees: "<<fileName1_<<"; "<<treeName1_
+ <<" and "<<fileName2_<<"; "<<treeName2_<<std::endl;
+ std::cout <<"Relative scaling factor histo1/histo2 set to "<<scaleFactor_<<std::endl;
+ std::cout <<"Minimum bin content is "<<minContent_<<std::endl;
+ std::cout <<"Number of free parameters is "<<nParams_<<std::endl;
+
+ g_data_smartptr = d1.Filter(Form("iExpt==%s",iExpt_.Data())).Graph(xName1_.Data(),yName1_.Data());
+ g_data_ = dynamic_cast<TGraph*>( g_data_smartptr.GetPtr() );
- Double_t* xs = new Double_t[nEntries];
- Double_t* ys = new Double_t[nEntries];
-
- for ( Int_t i=0; i < nEntries; ++i ) {
- tree1->GetEntry( i );
- xs[i] = x;
- ys[i] = y;
- }
+ // Get toy data
+ g_toy_smartptr = d2.Range(0,5000).Graph(xName1_.Data(),yName1_.Data());
+ g_toy_ = g_toy_smartptr.GetPtr();
theHisto_ = new TH2Poly("theHisto_", "", xMin_, xMax_, yMin_, yMax_);
//select the number of divisions to get us closest to minContent entries per bin
std::vector<Int_t> divisions;
- this->pickBinning(xs,ys,nEntries,divisions);
+ this->pickBinning(g_data_->GetX(),g_data_->GetY(),g_data_->GetN(),divisions);
//perform the adaptive bining based on histo1_
- this->getHisto(xMin_, xMax_, yMin_, yMax_, xs, ys, nEntries, divisions);
+ this->getHisto(xMin_, xMax_, yMin_, yMax_, g_data_->GetX(),g_data_->GetY(),g_data_->GetN(), divisions);
histo1_ = dynamic_cast<TH2Poly*>(theHisto_->Clone("histo1_"));
histo2_ = dynamic_cast<TH2Poly*>(theHisto_->Clone("histo2_"));
pullHisto_ = dynamic_cast<TH2Poly*>(theHisto_->Clone("pullHisto_"));
chiSqHisto_ = dynamic_cast<TH2Poly*>(theHisto_->Clone("chiSqHisto_"));
chiSqSignedHisto_ = dynamic_cast<TH2Poly*>(theHisto_->Clone("chiSqSignedHisto_"));
- delete[] xs;
- delete[] ys;
-
- //fill the two histograms from the trees
- TString drawString1, drawString2, weightString2;
- drawString1 += yName1_;
- drawString1 += ":";
- drawString1 += xName1_;
- drawString1 += ">>histo1_";
- drawString2 += yName2_;
- drawString2 += ":";
- drawString2 += xName2_;
- drawString2 += ">>histo2_";
- weightString2 += scaleFactor_;
-
- tree1->Draw(drawString1);
- tree2->Draw(drawString2,weightString2);
+
+ auto fillPoly = [this] ( const Double_t x, const Double_t y ) {
+ histo1_->Fill( x,y );
+ };
+ d1.Foreach( fillPoly, { xName1_.Data(), yName1_.Data()} );
+
+ ROOT::DisableImplicitMT(); // This is not thread-safe (thanks to TH2Poly)
+ auto fillPolyWeighted = [this] ( const Double_t x, const Double_t y ) {
+ histo2_->Fill( x,y,scaleFactor_ );
+ };
+ d2.Foreach( fillPolyWeighted, { xName1_.Data(), yName1_.Data()} );
histo1_->SetDirectory(0);
histo2_->SetDirectory(0);
pullHisto_->SetDirectory(0);
chiSqHisto_->SetDirectory(0);
chiSqSignedHisto_->SetDirectory(0);
-
- // close the file(s) containing the trees
- if (file1 != 0) {file1->Close();}
- delete file1;
- if (file2 != 0) {file2->Close();}
- delete file2;
}
void LauCalcChiSq::pickBinning(const Double_t* xs, const Double_t* ys, const Int_t nEntries, std::vector<Int_t>& divisions)

File Metadata

Mime Type
text/plain
Expires
Mon, Sep 29, 8:35 PM (22 h, 12 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
6546987
Default Alt Text
D55.1759174522.diff (11 KB)

Event Timeline

HEPForge.org · .