Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TMVA_SOFIE_GNN_Application.C
Go to the documentation of this file.
1// Macor evaluating a GNN model which was generated with the Parser macro
2//
3
4// need to add include path to find generated model file
5#ifdef __CLING__
7R__ADD_INCLUDE_PATH($ROOTSYS/runtutorials)
8#endif
9
10#include "encoder.hxx"
11#include "core.hxx"
12#include "decoder.hxx"
13#include "output_transform.hxx"
14
15#include "TMVA/SOFIE_common.hxx"
16#include "TRandom3.h"
17#include "TH1.h"
18#include "TCanvas.h"
19#include "TFile.h"
20#include "TTree.h"
21#include "TSystem.h"
22#include "ROOT/RDataFrame.hxx"
23
24using namespace TMVA::Experimental;
25using namespace TMVA::Experimental::SOFIE;
26
27double check_mem(std::string s = ""){
29 printf("%s - ",s.c_str());
31 printf(" Rmem = %8.3f MB, Vmem = %8.f3 MB \n",
32 p.fMemResident /1024., /// convert memory from kB to MB
33 p.fMemVirtual /1024.
34 );
35 return p.fMemResident / 1024.;
36}
37
38
39template<class T>
41 std::cout << " shape : " << ConvertShapeToString(t.GetShape()) << " size : " << t.GetSize() << "\n";
42 auto & shape = t.GetShape();
43 auto p = t.GetData();
44 size_t nrows = (shape.size() > 1) ? shape[0] : 1;
45 size_t ncols = (shape.size() > 1) ? t.GetStrides()[0] : shape[0];
46 for (size_t i = 0; i < nrows; i++) {
47 for (size_t j = 0; j < ncols; j++) {
48 if (j==ncols-1) {
49 if (j>10) std::cout << "... ";
50 std::cout << *p << std::endl;
51 }
52 else if (j<10)
53 std::cout << *p << ", ";
54 p++;
55 }
56 }
57 std::cout << std::endl;
58}
59void Print(GNN_Data & d, std::string txt = "") {
60 if (!txt.empty()) std::cout << std::endl << txt << std::endl;
61 std::cout << "node data:"; PrintTensor(d.node_data);
62 std::cout << "edge data:"; PrintTensor(d.edge_data);
63 std::cout << "global data:"; PrintTensor(d.global_data);
64 std::cout << "edge index:"; PrintTensor(d.edge_index);
65}
66
67struct SOFIE_GNN {
68 bool verbose = false;
69 TMVA_SOFIE_encoder::Session encoder;
70 TMVA_SOFIE_core::Session core;
71 TMVA_SOFIE_decoder::Session decoder;
72 TMVA_SOFIE_output_transform::Session output_transform;
73
74 std::vector<GNN_Data> Infer(const GNN_Data & data, int nsteps) {
75 // infer function
76 auto input_data = Copy(data);
77 if (verbose) Print(input_data,"input_data");
78 encoder.infer(input_data);
79 // latent0 is result of encoder. Need to copy because this stays the same
80 auto latent0 = Copy(input_data);
81 GNN_Data latent = input_data; // this can be a view
82 std::vector<GNN_Data> outputData;
83 for (int i = 0; i < nsteps; i++) {
84 if (verbose) Print(latent0,"input decoded data");
85 if (verbose) Print(latent,"latent data");
86 auto core_input = Concatenate(latent0, latent,1);
87 if (verbose) Print(core_input, "after concatenate");
88 core.infer(core_input);
89 if (verbose) Print(core_input, "after core inference");
90 // here I need to copy
91 latent = Copy(core_input);
92 decoder.infer(core_input);
93 output_transform.infer(core_input);
94 outputData.push_back(Copy(core_input));
95 }
96 return outputData;
97 }
98
99 SOFIE_GNN(bool v = false) : verbose(v) {}
100
101};
102
103const int num_max_nodes = 10;
104const int num_max_edges = 30;
105const int NODE_FEATURE_SIZE = 4;
106const int EDGE_FEATURE_SIZE = 4;
108
109std::vector<GNN_Data> GenerateData(int nevts, int seed) {
110 TRandom3 r(seed);
111 std::vector<GNN_Data> dataSet;
112 dataSet.reserve(nevts);
113 for (int i = 0; i < nevts; i++) {
114 // generate first number of nodes and edges
115 // size_t num_nodes = num_max_nodes;//r.Integer(num_max_nodes-2) + 2;
116 // size_t num_edges = num_max_edges; //r.Integer(num_max_edges-1) + 1;
117 size_t num_nodes = r.Integer(num_max_nodes-2) + 2;
118 size_t num_edges = r.Integer(num_max_edges-1) + 1;
119 GNN_Data gd;
120 gd.node_data = RTensor<float>({num_nodes, NODE_FEATURE_SIZE});
121 gd.edge_data = RTensor<float>({num_edges, EDGE_FEATURE_SIZE});
123 gd.edge_index = RTensor<int>({2, num_edges});
124
125 auto genValue = [&]() { return r.Rndm()*10 -5; };
126 auto genLink = [&] () { return r.Integer(num_nodes);};
127 std::generate(gd.node_data.begin(), gd.node_data.end(), genValue);
128 std::generate(gd.edge_data.begin(), gd.edge_data.end(), genValue);
129 std::generate(gd.global_data.begin(), gd.global_data.end(), genValue);
130 std::generate(gd.edge_index.begin(), gd.edge_index.end(), genLink);
131 dataSet.emplace_back(gd);
132 }
133 return dataSet;
134}
135
136std::vector<GNN_Data> ReadData(std::string treename, std::string filename) {
137 bool verbose = false;
138 ROOT::RDataFrame df(treename,filename);
139 auto ndata = df.Take<ROOT::RVec<float>>("node_data");
140 auto edata = df.Take<ROOT::RVec<float>>("edge_data");
141 auto gdata = df.Take<ROOT::RVec<float>>("global_data");
142 auto rdata = df.Take<ROOT::RVec<int>>("receivers");
143 auto sdata = df.Take<ROOT::RVec<int>>("senders");
144 int nevts = ndata.GetPtr()->size();
145 std::vector<GNN_Data> dataSet;
146 dataSet.reserve(nevts);
147 for (int i = 0; i < nevts; i++) {
148 GNN_Data gd;
149 auto & n = (*(ndata.GetPtr()))[i];
150 size_t num_nodes = n.size()/NODE_FEATURE_SIZE;
151 auto & e = (*(edata.GetPtr()))[i];
152 size_t num_edges = e.size()/EDGE_FEATURE_SIZE;
153 auto & g = (*(gdata.GetPtr()))[i];
154 gd.node_data = RTensor<float>(n.data(), {num_nodes, NODE_FEATURE_SIZE});
155 gd.edge_data = RTensor<float>(e.data(), {num_edges, EDGE_FEATURE_SIZE});
156 gd.global_data = RTensor<float>(g.data(), {1, GLOBAL_FEATURE_SIZE});
157 gd.edge_index = RTensor<int>({2, num_edges});
158 auto & r = (*(rdata.GetPtr()))[i];
159 auto & s = (*(sdata.GetPtr()))[i];
160 // need to copy receivers/senders in edge_index tensor
161 std::copy(r.begin(), r.end(), gd.edge_index.GetData());
162 std::copy(s.begin(), s.end(), gd.edge_index.GetData()+num_edges);
163
164 dataSet.emplace_back(Copy(gd)); // need to copy data in vector to own
165 if (i < 1 && verbose) Print(dataSet[i],"Input for Event" + std::to_string(i));
166 }
167 return dataSet;
168}
169
170
171void TMVA_SOFIE_GNN_Application (bool verbose = false)
172{
173 check_mem("Initial memory");
174 SOFIE_GNN gnn;
175 check_mem("After creating GNN");
176
177
178 const int seed = 111;
179 const int nproc_steps = 5;
180 // generate the input data
181
182 int nevts;
183 //std::cout << "generating data\n";
184 //nevts = 100;
185 //auto inputData = GenerateData(nevts, seed);
186
187 std::cout << "reading data\n";
188 auto inputData = ReadData("gdata","graph_data.root");
189 nevts = inputData.size();
190
191 //std::cout << "padding data\n";
192 //PadData(inputData) ;
193
194 auto h1 = new TH1D("h1","SOFIE Node data",40,1,0);
195 auto h2 = new TH1D("h2","SOFIE Edge data",40,1,0);
196 auto h3 = new TH1D("h3","SOFIE Global data",40,1,0);
197 std::cout << "doing inference...\n";
198
199
200 check_mem("Before evaluating");
201 TStopwatch w; w.Start();
202 for (int i = 0; i < nevts; i++) {
203 auto result = gnn.Infer(inputData[i], nproc_steps);
204 // compute resulting mean and plot them
205 auto & lr = result.back();
206 if (i < 1 && verbose) Print(lr,"Output for Event" + std::to_string(i));
207 h1->Fill(TMath::Mean(lr.node_data.begin(), lr.node_data.end()));
208 h2->Fill(TMath::Mean(lr.edge_data.begin(), lr.edge_data.end()));
209 h3->Fill(TMath::Mean(lr.global_data.begin(), lr.global_data.end()));
210 }
211 w.Stop();
212 w.Print();
213 check_mem("End evaluation");
214 auto c1 = new TCanvas("c1","SOFIE Results");
215 c1->Divide(1,3);
216 c1->cd(1); h1->Draw();
217 c1->cd(2); h2->Draw();
218 c1->cd(3); h3->Draw();
219
220
221 // compare with the reference
222 auto c2 = new TCanvas("c2","Reference Results");
223 auto file = TFile::Open("graph_data.root");
224 auto o1 = file->Get("h1");
225 auto o2 = file->Get("h2");
226 auto o3 = file->Get("h3");
227 c2->Divide(1,3);
228 c2->cd(1); o1->Draw();
229 c2->cd(2); o2->Draw();
230 c2->cd(3); o3->Draw();
231
232}
233
234int main () {
235
237}
238
239
#define d(i)
Definition RSha256.hxx:102
#define g(i)
Definition RSha256.hxx:105
#define e(i)
Definition RSha256.hxx:103
#define R__ADD_INCLUDE_PATH(PATH)
Definition Rtypes.h:497
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
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 filename
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 r
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 result
void Print(GNN_Data &d, std::string txt="")
const int EDGE_FEATURE_SIZE
const int GLOBAL_FEATURE_SIZE
const int num_max_edges
std::vector< GNN_Data > ReadData(std::string treename, std::string filename)
void PrintTensor(RTensor< T > &t)
const int num_max_nodes
void TMVA_SOFIE_GNN_Application(bool verbose=false)
const int NODE_FEATURE_SIZE
double check_mem(std::string s="")
std::vector< GNN_Data > GenerateData(int nevts, int seed)
R__EXTERN TSystem * gSystem
Definition TSystem.h:561
RResultPtr< COLL > Take(std::string_view column="")
Return a collection of values of a column (lazy action, returns a std::vector by default).
ROOT's RDataFrame offers a modern, high-level interface for analysis of data stored in TTree ,...
A "std::vector"-like collection of values implementing handy operation to analyse them.
Definition RVec.hxx:1529
The Canvas class.
Definition TCanvas.h:23
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:4089
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:670
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3344
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3066
RTensor is a container with contiguous memory and shape information.
Definition RTensor.hxx:162
Iterator end() noexcept
Definition RTensor.hxx:308
const Shape_t & GetStrides() const
Definition RTensor.hxx:243
std::size_t GetSize() const
Definition RTensor.hxx:241
Iterator begin() noexcept
Definition RTensor.hxx:305
const Shape_t & GetShape() const
Definition RTensor.hxx:242
Random number generator class based on M.
Definition TRandom3.h:27
Stopwatch class.
Definition TStopwatch.h:28
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
virtual int GetProcInfo(ProcInfo_t *info) const
Returns cpu and memory used by this process into the ProcInfo_t structure.
Definition TSystem.cxx:2489
return c1
Definition legend1.C:41
const Int_t n
Definition legend1.C:16
TH1F * h1
Definition legend1.C:5
return c2
Definition legend2.C:14
std::string ConvertShapeToString(std::vector< size_t > shape)
TMVA::Experimental::RTensor< T > Concatenate(TMVA::Experimental::RTensor< T > &t1, TMVA::Experimental::RTensor< T > &t2, int axis=0)
GNN_Data Copy(const GNN_Data &data)
Double_t Mean(Long64_t n, const T *a, const Double_t *w=nullptr)
Returns the weighted mean of an array a with length n.
Definition TMath.h:1089
TMVA_SOFIE_core::Session core
SOFIE_GNN(bool v=false)
TMVA_SOFIE_output_transform::Session output_transform
TMVA_SOFIE_decoder::Session decoder
TMVA_SOFIE_encoder::Session encoder
std::vector< GNN_Data > Infer(const GNN_Data &data, int nsteps)