 ROOT   Reference Guide
Searching...
No Matches
th2polyEurope.C File Reference

## Detailed Description  This tutorial illustrates how to create an histogram with polygonal bins (TH2Poly), fill it and draw it.

The initial data are stored in TMultiGraphs. They represent the european countries. The histogram filling is done according to a Mercator projection, therefore the bin contains should be proportional to the real surface of the countries.

The initial data have been downloaded from: http://www.maproom.psu.edu/dcw/ This database was developed in 1991/1992 and national boundaries reflect political reality as of that time.

The script is shooting npoints (script argument) randomly over the Europe area. The number of points inside the countries should be proportional to the country surface The estimated surface is compared to the surfaces taken from wikipedia.  Partitioning: Real Time = 1.32 seconds Cpu Time = 1.31 seconds
Filling : Real Time = 4.06 seconds Cpu Time = 4.06 seconds
THPoly Europe surface estimation error wrt wikipedia = 1.260096 per cent when using 500000 points
void th2polyEurope(Int_t npoints=500000)
{
Int_t i,j;
Double_t lon1 = -25;
Double_t lon2 = 35;
Double_t lat1 = 34;
Double_t lat2 = 72;
Double_t R = (lat2-lat1)/(lon2-lon1);
Int_t W = 800;
Int_t H = (Int_t)(R*800);
gStyle->SetStatX(0.28);
gStyle->SetStatY(0.45);
gStyle->SetStatW(0.15);
// Canvas used to draw TH2Poly (the map)
TCanvas *ce = new TCanvas("ce", "ce",0,0,W,H);
ce->SetGridx();
ce->SetGridy();
// Real surfaces taken from Wikipedia.
const Int_t nx = 36;
// see http://en.wikipedia.org/wiki/Area_and_population_of_European_countries
const char *countries[nx] = {
"france", "spain", "sweden", "germany", "finland",
"norway", "poland", "italy", "yugoslavia", "united_kingdom",
"romania", "belarus","greece", "czechoslovakia","bulgaria",
"iceland", "hungary","portugal","austria", "ireland",
"lithuania", "latvia", "estonia", "denmark", "netherlands",
"switzerland","moldova","belgium", "albania", "cyprus",
"luxembourg", "andorra","malta", "liechtenstein", "san_marino", "monaco" };
Float_t surfaces[nx] = {
547030, 505580, 449964, 357021, 338145,
324220, 312685, 301230, 255438, 244820,
237500, 207600, 131940, 127711, 110910,
103000, 93030, 89242, 83870, 70280,
65200, 64589, 45226, 43094, 41526,
41290, 33843, 30528, 28748, 9250,
2586, 468, 316, 160, 61, 2};
TH1F *h = new TH1F("h","Countries surfaces (in km^{2})",3,0,3);
for (i=0; i<nx; i++) h->Fill(countries[i], surfaces[i]);
h->LabelsDeflate();
TFile *f;
if (!f) {
printf("Cannot access europe.root. Is internet working ?\n");
return;
}
TH2Poly *p = new TH2Poly(
"Europe",
"Europe (bin contents are normalized to the surfaces in km^{2})",
lon1,lon2,lat1,lat2);
p->GetXaxis()->SetNdivisions(520);
p->GetXaxis()->SetTitle("longitude");
p->GetYaxis()->SetTitle("latitude");
p->SetContour(100);
TKey *key;
TIter nextkey(gDirectory->GetListOfKeys());
while ((key = (TKey*)nextkey())) {
if (obj->InheritsFrom("TMultiGraph")) {
mg = (TMultiGraph*)obj;
}
}
Double_t longitude, latitude;
Double_t x, y, pi4 = TMath::Pi()/4, alpha = TMath::Pi()/360;
gBenchmark->Start("Partitioning");
p->ChangePartition(100, 100);
gBenchmark->Show("Partitioning");
// Fill TH2Poly according to a Mercator projection.
gBenchmark->Start("Filling");
for (i=0; i<npoints; i++) {
longitude = r.Uniform(lon1,lon2);
latitude = r.Uniform(lat1,lat2);
x = longitude;
y = 38*TMath::Log(TMath::Tan(pi4+alpha*latitude));
p->Fill(x,y);
}
gBenchmark->Show("Filling");
Int_t nbins = p->GetNumberOfBins();
Double_t maximum = p->GetMaximum();
// h2 contains the surfaces computed from TH2Poly.
TH1F *h2 = (TH1F *)h->Clone("h2");
h2->Reset();
for (j=0; j<nx; j++) {
for (i=0; i<nbins; i++) {
if (strstr(countries[j],p->GetBinName(i+1))) {
h2->Fill(countries[j],p->GetBinContent(i+1));
h2->SetBinError(j, p->GetBinError(i+1));
}
}
}
// Normalize the TH2Poly bin contents to the real surfaces.
Double_t scale = surfaces/maximum;
for (i=0; i<nbins; i++) p->SetBinContent(i+1, scale*p->GetBinContent(i+1));
p->Draw("COLZ");
TCanvas *c1 = new TCanvas("c1", "c1",W+10,0,W-20,H);
c1->SetRightMargin(0.047);
scale = h->GetMaximum()/h2->GetMaximum();
h->SetStats(0);
h->SetLineColor(kRed-3);
h->SetLineWidth(2);
h->SetMarkerStyle(20);
h->SetMarkerColor(kBlue);
h->SetMarkerSize(0.8);
h->Draw("LP");
h->GetXaxis()->SetLabelFont(42);
h->GetXaxis()->SetLabelSize(0.03);
h->GetYaxis()->SetLabelFont(42);
h2->Scale(scale);
Double_t scale2=TMath::Sqrt(scale);
for (i=0; i<nx; i++) h2->SetBinError(i+1, scale2*h2->GetBinError(i+1));
h2->Draw("E SAME");
h2->SetMarkerStyle(20);
h2->SetMarkerSize(0.8);
TLegend *leg = new TLegend(0.5,0.67,0.92,0.8,nullptr,"NDC");
leg->SetTextFont(42);
leg->SetTextSize(0.025);
leg->AddEntry(h,"Real countries surfaces from Wikipedia (in km^{2})","lp");
leg->AddEntry(h2,"Countries surfaces from TH2Poly (with errors)","lp");
leg->Draw();
leg->Draw();
Double_t wikiSum = h->Integral();
Double_t polySum = h2->Integral();
Double_t error = TMath::Abs(wikiSum-polySum)/wikiSum;
printf("THPoly Europe surface estimation error wrt wikipedia = %f per cent when using %d points\n",100*error,npoints);
}
#define f(i)
Definition RSha256.hxx:104
#define h(i)
Definition RSha256.hxx:106
int Int_t
Definition RtypesCore.h:45
float Float_t
Definition RtypesCore.h:57
double Double_t
Definition RtypesCore.h:59
@ kRed
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
R__EXTERN TBenchmark * gBenchmark
Definition TBenchmark.h:59
#define gDirectory
Definition TDirectory.h:384
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
R__EXTERN TStyle * gStyle
Definition TStyle.h:433
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition TAttMarker.h:40
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition TAttMarker.h:45
virtual void Start(const char *name)
Starts Benchmark with the specified name.
virtual void Show(const char *name)
Stops Benchmark name and Prints results.
The Canvas class.
Definition TCanvas.h:23
virtual void ToggleEventStatus()
Toggle event statusbar.
Definition TCanvas.cxx:2425
A ROOT file is composed of a header, followed by consecutive data records (TKey instances) with a wel...
Definition TFile.h:53
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:4070
static Bool_t SetCacheFileDir(std::string_view cacheDir, Bool_t operateDisconnected=kTRUE, Bool_t forceCacheread=kFALSE)
Sets the directory where to locally stage/cache remote files.
Definition TFile.cxx:4606
1-D histogram with a float per channel (see TH1 documentation)}
Definition TH1.h:578
void Reset(Option_t *option="") override
Reset.
Definition TH1.cxx:9940
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.
Definition TH1.cxx:8890
virtual Double_t GetMaximum(Double_t maxval=FLT_MAX) const
Return maximum value smaller than maxval of bins in the range, unless the value has been overridden b...
Definition TH1.cxx:8372
virtual void SetBinError(Int_t bin, Double_t error)
Set the bin Error Note that this resets the bin eror option to be of Normal Type and for the non-empt...
Definition TH1.cxx:9033
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3325
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3047
virtual Double_t Integral(Option_t *option="") const
Return integral of bin contents.
Definition TH1.cxx:7812
virtual void Scale(Double_t c1=1, Option_t *option="")
Multiply this histogram by a constant c1.
Definition TH1.cxx:6553
TObject * Clone(const char *newname="") const override
Make a complete copy of the underlying object.
Definition TH1.cxx:2733
2D Histogram with Polygonal Bins
Definition TH2Poly.h:66
Book space in a file, create I/O buffers, to fill them, (un)compress them.
Definition TKey.h:28
To read a TObject* from the file.
Definition TKey.cxx:758
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition TMultiGraph.h:34
Mother of all ROOT objects.
Definition TObject.h:41
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition TObject.cxx:525
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:274
void SetGridx(Int_t value=1) override
void SetGridy(Int_t value=1) override
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition TStyle.cxx:1636
void SetStatX(Float_t x=0)
Definition TStyle.h:394
void SetStatW(Float_t w=0.19)
Definition TStyle.h:396
void SetStatY(Float_t y=0)
Definition TStyle.h:395
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
leg
Definition legend1.C:34
#define H(x, y, z)
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition TMath.h:756
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:662
constexpr Double_t Pi()
Definition TMath.h:37
Double_t Tan(Double_t)
Returns the tangent of an angle of x radians.
Definition TMath.h:600
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:123

Definition in file th2polyEurope.C.