Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
tprofile2polyRealistic.C File Reference

Detailed Description

View in nbviewer Open in SWAN
Different charges depending on region

#include <iostream>
#include <fstream>
#include <TProfile2Poly.h>
#include <TProfile2D.h>
#include <TCanvas.h>
#include <TRandom.h>
#include <TROOT.h>
void tprofile2polyRealistic(Int_t numEvents=100000)
{
int NUM_LS = 8;
TCanvas *c1 = new TCanvas("c1", "moving charge", 900, 400);
TCanvas *c2 = new TCanvas("c2", "Merge Individual moving charge plots", 800, 400);
c1->Divide(NUM_LS, 3);
c2->Divide(3,1);
// -------------------- Construct Reference plot bins ------------------------
auto new_avg = new TProfile2Poly();
auto tot_avg_ls = new TProfile2Poly[NUM_LS];
auto det_avg_ls = new TProfile2Poly[NUM_LS];
auto det_err_ls = new TProfile2Poly[NUM_LS];
auto tot_merge = new TProfile2Poly();
auto det_avg_merge = new TProfile2Poly();
auto det_err_merge = new TProfile2Poly();
float minx = -15;
float maxx = 15;
float miny = -15;
float maxy = 15;
float binsz = 0.5;
for (float i = minx; i < maxx; i += binsz) {
for (float j = miny; j < maxy; j += binsz) {
tot_merge->AddBin(i, j, i + binsz, j + binsz);
for (int l = 0; l < NUM_LS; ++l) {
tot_avg_ls[l].AddBin(i, j, i + binsz, j + binsz);
}
}
}
// -------------------- Construct detector bins ------------------------
auto h2p = new TH2Poly();
auto tp2p = new TProfile2Poly();
std::ifstream infile;
TString dir = gROOT->GetTutorialDir();
dir.Append("/hist/data/tprofile2poly_tutorial.data");
infile.open(dir.Data());
if (!infile) // Verify that the file was open successfully
{
std::cerr << dir.Data() << std::endl; // Report error
std::cerr << "Error code: " << std::strerror(errno) << std::endl; // Get some info as to why
return;
}
std::cout << " WE ARE AFTER LOADING DATA " << std::endl;
std::vector<std::pair<Double_t, Double_t>> allCoords;
while (infile >> a >> b)
allCoords.emplace_back(a, b);
if (allCoords.size() % 3 != 0) {
cout << "[ERROR] Bad file" << endl;
return;
}
Double_t x[3], y[3];
for (int i = 0; i < allCoords.size(); i += 3) {
x[0] = allCoords[i + 0].first;
y[0] = allCoords[i + 0].second;
x[1] = allCoords[i + 1].first;
y[1] = allCoords[i + 1].second;
x[2] = allCoords[i + 2].first;
y[2] = allCoords[i + 2].second;
det_avg_merge->AddBin(3, x, y);
det_err_merge->AddBin(3, x, y);
for (int l = 0; l < NUM_LS; ++l) {
det_avg_ls[l].AddBin(3, x, y);
det_err_ls[l].AddBin(3, x, y);
}
}
std::cout << " WE ARE AFTER ADDING BINS " << std::endl;
// -------------------- Simulate particles ------------------------
TRandom ran;
// moving error
Double_t xoffset1 = 0;
Double_t yoffset1 = 0;
Double_t xoffset2 = 0;
Double_t yoffset2 = 0;
for (int i = 0; i <= NUM_LS-1; ++i) { // LumiSection
std::cout << "[In Progress] LumiSection " << i << std::endl;
for (int j = 0; j < numEvents; ++j) { // Events
Double_t r1 = ran.Gaus(0, 10);
Double_t r2 = ran.Gaus(0, 8);
Double_t rok = ran.Gaus(10, 1);
Double_t rbad1 = ran.Gaus(8, 5);
Double_t rbad2 = ran.Gaus(-8, 5);
Double_t val = rok;
xoffset1 += 0.00002;
yoffset1 += 0.00002;
xoffset2 += 0.00003;
yoffset2 += 0.00004;
if (r2 > 3. - yoffset1 && r2 < 8. - yoffset1 &&
r1 > 1. + xoffset1 && r1 < 5. + xoffset1 ) {
val -= rbad1;
}
if (r2 > -10 + yoffset2 && r2 < -8 + yoffset2 &&
r1 > -6 + xoffset2 && r1 < 8 + xoffset2 ) {
val -= rbad2;
}
tot_avg_ls[i].Fill(r1, r2, val);
det_avg_ls[i].Fill(r1, r2, val);
det_err_ls[i].Fill(r1, r2, val);
}
std::string title;
c1->cd(i+1);
title = "Global View: Avg in LS " + std::to_string(i);
tot_avg_ls[i].SetTitle(title.c_str());
tot_avg_ls[i].SetStats(false);
tot_avg_ls[i].Draw("COLZ");
c1->Update();
c1->cd((i+1)+NUM_LS);
title = "Detector View: Avg in LS " + std::to_string(i);
det_avg_ls[i].SetTitle(title.c_str());
det_avg_ls[i].SetStats(false);
det_avg_ls[i].Draw("COLZ");
c1->Update();
c1->cd((i+1)+(NUM_LS*2));
title = "Detector View: Error in LS " + std::to_string(i);
det_err_ls[i].SetTitle(title.c_str());
det_err_ls[i].SetStats(false);
det_err_ls[i].SetContentToError();
det_err_ls[i].Draw("COLZ");
c1->Update();
}
std::vector<TProfile2Poly*> tot_avg_v;
std::vector<TProfile2Poly*> det_avg_v;
for (int t = 0; t < NUM_LS; t++){
tot_avg_v.push_back(&tot_avg_ls[t]);
det_avg_v.push_back(&det_avg_ls[t]);
}
std::cout << "[In Progress] Merging" << std::endl;
tot_merge->Merge(tot_avg_v);
c2->cd(1);
tot_merge->SetTitle("Total average merge");
tot_merge->Draw("COLZ");
det_avg_merge->Merge(det_avg_v);
c2->cd(2);
det_avg_merge->SetTitle("Detector average merge");
det_avg_merge->SetContentToAverage(); // implicit
det_avg_merge->Draw("COLZ");
det_err_merge->Merge(det_avg_v);
c2->cd(3);
det_err_merge->SetTitle("Detector error merge");
det_err_merge->SetContentToError();
det_err_merge->Draw("COLZ");
}
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
#define gROOT
Definition TROOT.h:406
The Canvas class.
Definition TCanvas.h:23
2D Histogram with Polygonal Bins
Definition TH2Poly.h:66
2D Profile Histogram with Polygonal Bins.
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition TRandom.cxx:275
Basic string class.
Definition TString.h:139
const char * Data() const
Definition TString.h:376
TString & Append(const char *cs)
Definition TString.h:572
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
return c2
Definition legend2.C:14
TLine l
Definition textangle.C:4
Author
Filip Ilic

Definition in file tprofile2polyRealistic.C.