Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RFieldProvider.hxx
Go to the documentation of this file.
1/*************************************************************************
2 * Copyright (C) 1995-2021, Rene Brun and Fons Rademakers. *
3 * All rights reserved. *
4 * *
5 * For the licensing terms see $ROOTSYS/LICENSE. *
6 * For the list of contributors see $ROOTSYS/README/CREDITS. *
7 *************************************************************************/
8
9#ifndef ROOT_Browsable_RFieldProvider
10#define ROOT_Browsable_RFieldProvider
11
12#include "TH1.h"
13#include "TMath.h"
14#include <map>
15#include <memory>
16#include <string>
17#include <utility>
18
20
22#include <ROOT/RPageStorage.hxx>
23#include <ROOT/RNTupleView.hxx>
24
25#include "RFieldHolder.hxx"
26
27using namespace ROOT::Browsable;
28
29using namespace std::string_literals;
30
31// FIXME: this exposes RField and RIntegralField into the global namespace
32template<typename T>
34template<typename T>
36
37// ==============================================================================================
38
39/** \class RFieldProvider
40\ingroup rbrowser
41\brief Base class for provider of RNTuple drawing
42\author Sergey Linev <S.Linev@gsi.de>
43\date 2021-03-09
44\warning This is part of the ROOT 7 prototype! It will change without notice. It might trigger earthquakes. Feedback is welcome!
45*/
46
47class RFieldProvider : public RProvider {
49 private:
50 std::shared_ptr<ROOT::Experimental::RNTupleReader> fNtplReader;
51 std::unique_ptr<TH1> fHist;
52
53 /** Test collected entries if it looks like integer values and one can use better binning */
55 {
56 auto len = fHist->GetBufferLength();
57 auto buf = fHist->GetBuffer();
58
59 if (!buf || (len < 5))
60 return;
61
62 Double_t min = buf[1], max = buf[1];
63 Bool_t is_integer = kTRUE;
64
65 for (Int_t n = 0; n < len; ++n) {
66 Double_t v = buf[2 + 2*n];
67 if (v > max) max = v;
68 if (v < min) min = v;
69 if (TMath::Abs(v - TMath::Nint(v)) > 1e-5) { is_integer = kFALSE; break; }
70 }
71
72 // special case when only integer values in short range - better binning
73 if (is_integer && (max-min < 100)) {
74 max += 2;
75 if (min > 1) min -= 2;
76 int npoints = TMath::Nint(max - min);
77 std::unique_ptr<TH1> h1 = std::make_unique<TH1F>(fHist->GetName(), fHist->GetTitle(), npoints, min, max);
78 h1->SetDirectory(nullptr);
79 for (Int_t n = 0; n < len; ++n)
80 h1->Fill(buf[2 + 2*n], buf[1 + 2*n]);
81 std::swap(fHist, h1);
82 }
83 }
84
85 template <typename T>
87 {
88 std::string title = "Drawing of RField "s + field.GetFieldName();
89
90 fHist = std::make_unique<TH1F>("hdraw", title.c_str(), 100, 0, 0);
91 fHist->SetDirectory(nullptr);
92
93 auto bufsize = (fHist->GetBufferSize() - 1) / 2;
94 int cnt = 0;
95 if (bufsize > 10) bufsize-=3; else bufsize = -1;
96
97 for (auto i : view.GetFieldRange()) {
98 fHist->Fill(view(i));
99 if (++cnt == bufsize) {
101 ++cnt;
102 }
103 }
104 if (cnt < bufsize)
106
107 fHist->BufferEmpty();
108 }
109
110 template<typename T>
112 {
113 auto view = fNtplReader->GetView<T>(field.GetOnDiskId());
114 FillHistogramImpl(field, view);
115 }
116
117 template<typename T>
118 void FillHistogram(const RField<T> &field)
119 {
120 auto view = fNtplReader->GetView<T>(field.GetOnDiskId());
121 FillHistogramImpl(field, view);
122 }
123
125 {
126 std::map<std::string, int> values;
127
128 int nentries = 0;
129
130 auto view = fNtplReader->GetView<std::string>(field.GetOnDiskId());
131 for (auto i : view.GetFieldRange()) {
132 std::string v = view(i);
133 nentries++;
134 auto iter = values.find(v);
135 if (iter != values.end())
136 iter->second++;
137 else if (values.size() >= 50)
138 return;
139 else
140 values[v] = 0;
141 }
142
143 // now create histogram with labels
144
145 std::string title = "Drawing of RField "s + field.GetFieldName();
146 fHist = std::make_unique<TH1F>("h",title.c_str(),3,0,3);
147 fHist->SetDirectory(nullptr);
148 fHist->SetStats(0);
149 fHist->SetEntries(nentries);
150 fHist->SetCanExtend(TH1::kAllAxes);
151 for (auto &entry : values)
152 fHist->Fill(entry.first.c_str(), entry.second);
153 fHist->LabelsDeflate();
154 fHist->Sumw2(kFALSE);
155 }
156
157 public:
158 explicit RDrawVisitor(std::shared_ptr<ROOT::Experimental::RNTupleReader> ntplReader) : fNtplReader(ntplReader) {}
159
161 return fHist.release();
162 }
163
164 void VisitField(const ROOT::Experimental::RFieldBase & /* field */) final {}
165 void VisitBoolField(const RField<bool> &field) final { FillHistogram(field); }
166 void VisitFloatField(const RField<float> &field) final { FillHistogram(field); }
167 void VisitDoubleField(const RField<double> &field) final { FillHistogram(field); }
168 void VisitCharField(const RField<char> &field) final { FillHistogram(field); }
169 void VisitInt8Field(const RIntegralField<std::int8_t> &field) final { FillHistogram(field); }
170 void VisitInt16Field(const RIntegralField<std::int16_t> &field) final { FillHistogram(field); }
171 void VisitInt32Field(const RIntegralField<std::int32_t> &field) final { FillHistogram(field); }
172 void VisitInt64Field(const RIntegralField<std::int64_t> &field) final { FillHistogram(field); }
173 void VisitStringField(const RField<std::string> &field) final { FillStringHistogram(field); }
177 void VisitUInt8Field(const RIntegralField<std::uint8_t> &field) final { FillHistogram(field); }
179 {
180 if (const auto f32 = field.As32Bit()) {
181 FillHistogram(*f32);
182 } else if (const auto f64 = field.As64Bit()) {
183 FillHistogram(*f64);
184 }
185 }
186 }; // class RDrawVisitor
187
188public:
189 // virtual ~RFieldProvider() = default;
190
192 {
193 if (!holder) return nullptr;
194
195 auto ntplReader = holder->GetNtplReader();
196 std::string name = holder->GetParentName();
197
198 const auto fieldName = ntplReader->GetDescriptor().GetFieldDescriptor(holder->GetId()).GetFieldName();
199 const auto qualifiedFieldName = ntplReader->GetDescriptor().GetQualifiedFieldName(holder->GetId());
200 auto view = ntplReader->GetView<void>(qualifiedFieldName);
201 name.append(fieldName);
202
203 RDrawVisitor drawVisitor(ntplReader);
204 view.GetField().AcceptVisitor(drawVisitor);
205 return drawVisitor.MoveHist();
206 }
207};
208
209#endif
#define e(i)
Definition RSha256.hxx:103
constexpr Bool_t kFALSE
Definition RtypesCore.h:94
constexpr Bool_t kTRUE
Definition RtypesCore.h:93
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 Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t UChar_t len
char name[80]
Definition TGX11.cxx:110
int nentries
auto GetParentName() const
auto GetId() const
auto GetNtplReader() const
void VisitCardinalityField(const ROOT::Experimental::RCardinalityField &field) final
void VisitUInt64Field(const RIntegralField< std::uint64_t > &field) final
void VisitInt8Field(const RIntegralField< std::int8_t > &field) final
void FillStringHistogram(const RField< std::string > &field)
std::unique_ptr< TH1 > fHist
void TestHistBuffer()
Test collected entries if it looks like integer values and one can use better binning.
void FillHistogram(const RField< T > &field)
void VisitUInt8Field(const RIntegralField< std::uint8_t > &field) final
void VisitDoubleField(const RField< double > &field) final
void VisitInt64Field(const RIntegralField< std::int64_t > &field) final
void VisitCharField(const RField< char > &field) final
void VisitUInt32Field(const RIntegralField< std::uint32_t > &field) final
void FillHistogramImpl(const ROOT::Experimental::RFieldBase &field, ROOT::Experimental::RNTupleView< T > &view)
void VisitField(const ROOT::Experimental::RFieldBase &) final
void VisitInt32Field(const RIntegralField< std::int32_t > &field) final
void FillHistogram(const RIntegralField< T > &field)
void VisitInt16Field(const RIntegralField< std::int16_t > &field) final
void VisitUInt16Field(const RIntegralField< std::uint16_t > &field) final
std::shared_ptr< ROOT::Experimental::RNTupleReader > fNtplReader
void VisitStringField(const RField< std::string > &field) final
void VisitBoolField(const RField< bool > &field) final
RDrawVisitor(std::shared_ptr< ROOT::Experimental::RNTupleReader > ntplReader)
void VisitFloatField(const RField< float > &field) final
Base class for provider of RNTuple drawing.
TH1 * DrawField(RFieldHolder *holder)
Provider of different browsing methods for supported classes.
Definition RProvider.hxx:37
Abstract base class for classes implementing the visitor design pattern.
An artificial field that transforms an RNTuple column that contains the offset of collections into co...
Definition RField.hxx:268
A field translates read and write calls from/to underlying columns to/from tree values.
const std::string & GetFieldName() const
DescriptorId_t GetOnDiskId() const
Classes with dictionaries that can be inspected by TClass.
Definition RField.hxx:241
RNTupleGlobalRange GetFieldRange() const
An RNTupleView for a known type.
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:59
virtual void SetDirectory(TDirectory *dir)
By default, when a histogram is created, it is added to the list of histogram objects in the current ...
Definition TH1.cxx:8970
@ kAllAxes
Definition TH1.h:76
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3346
const Int_t n
Definition legend1.C:16
TH1F * h1
Definition legend1.C:5
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
Definition TMath.h:697
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:123