Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
ntpl009_parallelWriter.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_ntuple
3/// \notebook
4/// Example of multi-threaded writes using RNTupleParallelWriter. Adapted from the ntpl007_mtFill tutorial.
5///
6/// \macro_image
7/// \macro_code
8///
9/// \date Feburary 2024
10/// \author The ROOT Team
11
12// NOTE: The RNTupleParallelWriter class is experimental at this point.
13// Functionality and interface are still subject to changes.
14
16#include <ROOT/RNTupleModel.hxx>
19
20#include <TCanvas.h>
21#include <TH1F.h>
22#include <TH2F.h>
23#include <TRandom.h>
24#include <TRandom3.h>
25#include <TStyle.h>
26#include <TSystem.h>
27
28#include <atomic>
29#include <memory>
30#include <mutex>
31#include <thread>
32#include <vector>
33#include <utility>
34
35// Import classes from experimental namespace for the time being
37
38// Where to store the ntuple of this example
39constexpr char const *kNTupleFileName = "ntpl009_parallelWriter.root";
40
41// Number of parallel threads to fill the ntuple
42constexpr int kNWriterThreads = 4;
43
44// Number of events to generate is kNEventsPerThread * kNWriterThreads
45constexpr int kNEventsPerThread = 25000;
46
47// Thread function to generate and write events
48void FillData(RNTupleParallelWriter *writer)
49{
50 static std::atomic<std::uint32_t> gThreadId;
51 const auto threadId = ++gThreadId;
52
53 auto prng = std::make_unique<TRandom3>();
54 prng->SetSeed();
55
56 auto fillContext = writer->CreateFillContext();
57 auto entry = fillContext->CreateEntry();
58
59 auto id = entry->GetPtr<std::uint32_t>("id");
60 *id = threadId;
61 auto vpx = entry->GetPtr<std::vector<float>>("vpx");
62 auto vpy = entry->GetPtr<std::vector<float>>("vpy");
63 auto vpz = entry->GetPtr<std::vector<float>>("vpz");
64
65 for (int i = 0; i < kNEventsPerThread; i++) {
66 vpx->clear();
67 vpy->clear();
68 vpz->clear();
69
70 int npx = static_cast<int>(prng->Rndm(1) * 15);
71 // Set the field data for the current event
72 for (int j = 0; j < npx; ++j) {
73 float px, py, pz;
74 prng->Rannor(px, py);
75 pz = px * px + py * py;
76
77 vpx->emplace_back(px);
78 vpy->emplace_back(py);
79 vpz->emplace_back(pz);
80 }
81
82 fillContext->Fill(*entry);
83 }
84}
85
86// Generate kNEvents with multiple threads in kNTupleFileName
87void Write()
88{
89 // Create the data model
90 auto model = ROOT::RNTupleModel::CreateBare();
91 model->MakeField<std::uint32_t>("id");
92 model->MakeField<std::vector<float>>("vpx");
93 model->MakeField<std::vector<float>>("vpy");
94 model->MakeField<std::vector<float>>("vpz");
95
96 // Create RNTupleWriteOptions to make the writing commit multiple clusters (so that "Entry Id vs Thread Id" shows the
97 // interleaved clusters).
99 options.SetApproxZippedClusterSize(1024 * 1024);
100
101 // We hand-over the data model to a newly created ntuple of name "NTuple", stored in kNTupleFileName
102 auto writer = RNTupleParallelWriter::Recreate(std::move(model), "NTuple", kNTupleFileName, options);
103
104 std::vector<std::thread> threads;
105 for (int i = 0; i < kNWriterThreads; ++i)
106 threads.emplace_back(FillData, writer.get());
107 for (int i = 0; i < kNWriterThreads; ++i)
108 threads[i].join();
109
110 // The writer unique pointer goes out of scope here. On destruction, the writer flushes unwritten data to disk
111 // and closes the attached ROOT file.
112}
113
114// For all of the events, histogram only one of the written vectors
115void Read()
116{
118 auto viewVpx = reader->GetView<float>("vpx._0");
119
120 gStyle->SetOptStat(0);
121
122 TCanvas *c1 = new TCanvas("c2", "Multi-Threaded Filling Example", 200, 10, 1500, 500);
123 c1->Divide(2, 1);
124
125 c1->cd(1);
126 TH1F h("h", "This is the px distribution", 100, -4, 4);
127 h.SetFillColor(48);
128 // Iterate through all values of vpx in all events
129 for (auto i : viewVpx.GetFieldRange())
130 h.Fill(viewVpx(i));
131 // Prevent the histogram from disappearing
132 h.DrawCopy();
133
134 c1->cd(2);
135 auto nEvents = reader->GetNEntries();
136 auto viewId = reader->GetView<std::uint32_t>("id");
137 TH2F hFillSequence("", "Entry Id vs Thread Id;Entry Sequence Number;Filling Thread", 100, 0, nEvents, 100, 0,
138 kNWriterThreads + 1);
139 for (auto i : reader->GetEntryRange())
140 hFillSequence.Fill(i, viewId(i));
141 hFillSequence.DrawCopy();
142}
143
145{
146 Write();
147 Read();
148}
#define h(i)
Definition RSha256.hxx:106
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
R__EXTERN TStyle * gStyle
Definition TStyle.h:442
A writer to fill an RNTuple from multiple contexts.
static std::unique_ptr< RNTupleModel > CreateBare()
Creates a "bare model", i.e. an RNTupleModel with no default entry.
static std::unique_ptr< RNTupleReader > Open(std::string_view ntupleName, std::string_view storage, const ROOT::RNTupleReadOptions &options=ROOT::RNTupleReadOptions())
Open an RNTuple for reading.
Common user-tunable settings for storing RNTuples.
The Canvas class.
Definition TCanvas.h:23
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:877
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:307
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:1642
return c1
Definition legend1.C:41
void FillData(BinData &dv, const TH1 *hist, TF1 *func=nullptr)
fill the data vector from a TH1.
ROOT::RNTupleGlobalRange GetFieldRange(const ROOT::RFieldBase &field, const ROOT::Internal::RPageSource &pageSource)
Helper to get the iteration space of the given field that needs to be connected to the given page sou...