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

Detailed Description

View in nbviewer Open in SWAN
Work with befriended RNTuples.

Adapted from tree3.C

// NOTE: The RNTuple classes are experimental at this point.
// Functionality, interface, and data format is still subject to changes.
// Do not use for real data!
#include <TCanvas.h>
#include <TH1F.h>
#include <TMath.h>
#include <TRandom.h>
#include <vector>
constexpr char const* kNTupleMainFileName = "ntpl006_data.root";
constexpr char const* kNTupleFriendFileName = "ntpl006_reco.root";
using RNTupleModel = ROOT::Experimental::RNTupleModel;
using RNTupleReader = ROOT::Experimental::RNTupleReader;
using RNTupleWriter = ROOT::Experimental::RNTupleWriter;
void Generate()
{
constexpr int kMaxTrack = 500;
auto modelData = RNTupleModel::Create();
auto fldPx = modelData->MakeField<std::vector<float>>("px");
auto fldPy = modelData->MakeField<std::vector<float>>("py");
auto fldPz = modelData->MakeField<std::vector<float>>("pz");
auto modelFriend = RNTupleModel::Create();
auto fldPt = modelFriend->MakeField<std::vector<float>>("pt");
auto ntupleData = RNTupleWriter::Recreate(std::move(modelData), "data", kNTupleMainFileName);
auto ntupleReco = RNTupleWriter::Recreate(std::move(modelFriend), "reco", kNTupleFriendFileName);
for (int i=0; i < 1000; i++) {
int ntracks = gRandom->Rndm() * (kMaxTrack - 1);
for (int n = 0; n < ntracks; n++) {
fldPx->emplace_back(gRandom->Gaus( 0, 1));
fldPy->emplace_back(gRandom->Gaus( 0, 2));
fldPz->emplace_back(gRandom->Gaus(10, 5));
fldPt->emplace_back(TMath::Sqrt(fldPx->at(n) * fldPx->at(n) + fldPy->at(n) * fldPy->at(n)));
}
ntupleData->Fill();
ntupleReco->Fill();
fldPx->clear();
fldPy->clear();
fldPz->clear();
fldPt->clear();
}
}
void ntpl006_friends()
{
Generate();
std::vector<RNTupleReader::ROpenSpec> friends{
{"data", kNTupleMainFileName},
{"reco", kNTupleFriendFileName} };
auto ntuple = RNTupleReader::OpenFriends(friends);
auto c = new TCanvas("c", "", 200, 10, 700, 500);
TH1F h("h", "pz {pt > 3.}", 100, -15, 35);
auto viewPz = ntuple->GetView<float>("data.pz._0");
auto viewPt = ntuple->GetView<float>("reco.pt._0");
for (auto i : viewPt.GetFieldRange()) {
if (viewPt(i) > 3.)
h.Fill(viewPz(i));
}
h.SetFillColor(48);
h.DrawCopy();
}
#define c(i)
Definition RSha256.hxx:101
#define h(i)
Definition RSha256.hxx:106
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
The RNTupleModel encapulates the schema of an ntuple.
An RNTuple that is used to read data from storage.
An RNTuple that gets filled with entries (data) and writes them to storage.
The Canvas class.
Definition TCanvas.h:23
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:621
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
Double_t Rndm() override
Machine independent random number generator.
Definition TRandom.cxx:559
const Int_t n
Definition legend1.C:16
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:662
Date
April 2020
Author
The ROOT Team

Definition in file ntpl006_friends.C.