Write and read STL vectors with RNTuple. Adapted from the hvector tree tutorial.
Author: The ROOT Team
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, April 17, 2024 at 11:24 AM.
%%cpp -d
#include <ROOT/RNTupleModel.hxx>
#include <ROOT/RNTupleReader.hxx>
#include <ROOT/RNTupleWriter.hxx>
#include <TCanvas.h>
#include <TH1F.h>
#include <TRandom.h>
#include <TSystem.h>
#include <cstdio>
#include <iostream>
#include <memory>
#include <vector>
#include <utility>
using RNTupleModel = ROOT::Experimental::RNTupleModel;
using RNTupleReader = ROOT::Experimental::RNTupleReader;
using RNTupleWriter = ROOT::Experimental::RNTupleWriter;
constexpr char const* kNTupleFileName = "ntpl002_vector.root";
constexpr int kUpdateGuiFreq = 1000;
constexpr int kNEvents = 25000;
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! Import classes from experimental namespace for the time being Where to store the ntuple of this example Update the histogram GUI every so many fills Number of events to generate Generate kNEvents with vectors in kNTupleFileName
%%cpp -d
void Write()
{
// We create a unique pointer to an empty data model
auto model = RNTupleModel::Create();
// Creating fields of std::vector is the same as creating fields of simple types. As a result, we get
// shared pointers of the given type
std::shared_ptr<std::vector<float>> fldVpx = model->MakeField<std::vector<float>>("vpx");
auto fldVpy = model->MakeField<std::vector<float>>("vpy");
auto fldVpz = model->MakeField<std::vector<float>>("vpz");
auto fldVrand = model->MakeField<std::vector<float>>("vrand");
// We hand-over the data model to a newly created ntuple of name "F", stored in kNTupleFileName
// In return, we get a unique pointer to an ntuple that we can fill
auto ntuple = RNTupleWriter::Recreate(std::move(model), "F", kNTupleFileName);
TH1F hpx("hpx", "This is the px distribution", 100, -4, 4);
hpx.SetFillColor(48);
auto c1 = new TCanvas("c1", "Dynamic Filling Example", 200, 10, 700, 500);
gRandom->SetSeed();
for (int i = 0; i < kNEvents; i++) {
int npx = static_cast<int>(gRandom->Rndm(1) * 15);
fldVpx->clear();
fldVpy->clear();
fldVpz->clear();
fldVrand->clear();
// Set the field data for the current event
for (int j = 0; j < npx; ++j) {
float px, py, pz;
gRandom->Rannor(px, py);
pz = px*px + py*py;
auto random = gRandom->Rndm(1);
hpx.Fill(px);
fldVpx->emplace_back(px);
fldVpy->emplace_back(py);
fldVpz->emplace_back(pz);
fldVrand->emplace_back(random);
}
// Gui updates
if (i && (i % kUpdateGuiFreq) == 0) {
if (i == kUpdateGuiFreq) hpx.Draw();
c1->Modified();
c1->Update();
if (gSystem->ProcessEvents())
break;
}
ntuple->Fill();
}
hpx.DrawCopy();
// The ntuple unique pointer goes out of scope here. On destruction, the ntuple flushes unwritten data to disk
// and closes the attached ROOT file.
}
For all of the events, histogram only one of the written vectors
%%cpp -d
void Read()
{
// Get a unique pointer to an empty RNTuple model
auto model = RNTupleModel::Create();
// We only define the fields that are needed for reading
auto fldVpx = model->MakeField<std::vector<float>>("vpx");
// Create an ntuple without imposing a specific data model. We could generate the data model from the ntuple
// but here we prefer the view because we only want to access a single field
auto ntuple = RNTupleReader::Open(std::move(model), "F", kNTupleFileName);
// Quick overview of the ntuple's key meta-data
ntuple->PrintInfo();
std::cout << "Entry number 42 in JSON format:" << std::endl;
ntuple->Show(41);
// In a future version of RNTuple, there will be support for ntuple->Scan()
TCanvas *c2 = new TCanvas("c2", "Dynamic Filling Example", 200, 10, 700, 500);
TH1F h("h", "This is the px distribution", 100, -4, 4);
h.SetFillColor(48);
// Iterate through all the events using i as event number and as an index for accessing the view
for (auto entryId : *ntuple) {
ntuple->LoadEntry(entryId);
for (auto px : *fldVpx) {
h.Fill(px);
}
if (entryId && (entryId % kUpdateGuiFreq) == 0) {
if (entryId == kUpdateGuiFreq) h.Draw();
c2->Modified();
c2->Update();
if (gSystem->ProcessEvents())
break;
}
}
// Prevent the histogram from disappearing
h.DrawCopy();
}
Write();
Read();
************************************ NTUPLE ************************************ * N-Tuple : F * * Entries : 25000 * ******************************************************************************** * Field 1 : vpx (std::vector<float>) * * Field 1.1 : _0 (float) * * Field 2 : vpy (std::vector<float>) * * Field 2.1 : _0 (float) * * Field 3 : vpz (std::vector<float>) * * Field 3.1 : _0 (float) * * Field 4 : vrand (std::vector<float>) * * Field 4.1 : _0 (float) * ******************************************************************************** Entry number 42 in JSON format: { "vpx": [-1.19172, 0.904947, -1.07199, 0.33103, 1.63725, -0.867451], "vpy": [-0.418819, 0.972805, -1.29367, 1.37038, 1.00845, -0.920446], "vpz": [1.59561, 1.76528, 2.82273, 1.98753, 3.69755, 1.59969], "vrand": [0.32247, 0.584583, 0.594393, 0.844586, 0.221563, 0.146908] }
Warning in <[ROOT.NTuple] Warning /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master/tree/ntuple/v7/src/RPageStorageFile.cxx:56 in ROOT::Experimental::Internal::RPageSinkFile::RPageSinkFile(std::string_view, const ROOT::Experimental::RNTupleWriteOptions&)::<lambda()>>: The RNTuple file format will change. Do not store real data with this version of RNTuple! Warning in <[ROOT.NTuple] Warning /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master/tree/ntuple/v7/src/RPageStorageFile.cxx:262 in ROOT::Experimental::Internal::RPageSourceFile::InitDescriptor(const ROOT::Experimental::RNTuple&)::<lambda()>>: Pre-release format version: RC 2
Draw all canvases
gROOT->GetListOfCanvases()->Draw()