Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RSofieReader.hxx
Go to the documentation of this file.
1/**********************************************************************************
2 * Project: ROOT - a Root-integrated toolkit for multivariate data analysis *
3 * Package: TMVA * *
4 * *
5 * Description: *
6 * *
7 * Authors: *
8 * Lorenzo Moneta *
9 * *
10 * Copyright (c) 2022: *
11 * CERN, Switzerland *
12 * *
13 **********************************************************************************/
14
15
16#ifndef TMVA_RSOFIEREADER
17#define TMVA_RSOFIEREADER
18
19
20#include <string>
21#include <vector>
22#include <memory> // std::unique_ptr
23#include <sstream> // std::stringstream
24#include <iostream>
25#include "TROOT.h"
26#include "TSystem.h"
27#include "TError.h"
28#include "TInterpreter.h"
29#include "TUUID.h"
30#include "TMVA/RTensor.hxx"
31#include "Math/Util.h"
32
33namespace TMVA {
34namespace Experimental {
35
36
37
38
39/// TMVA::RSofieReader class for reading external Machine Learning models
40/// in ONNX files, Keras .h5 files or PyTorch .pt files
41/// and performing the inference using SOFIE
42/// It is reccomended to use ONNX if possible since there is a larger support for
43/// model operators.
44
46
47
48public:
49 /// Dummy constructor which needs model loading afterwards
51 /// Create TMVA model from ONNX file
52 /// print level can be 0 (minimal) 1 with info , 2 with all ONNX parsing info
53 RSofieReader(const std::string &path, std::vector<std::vector<size_t>> inputShapes = {}, int verbose = 0)
54 {
55 Load(path, inputShapes, verbose);
56 }
57
58 void Load(const std::string &path, std::vector<std::vector<size_t>> inputShapes = {}, int verbose = 0)
59 {
60
61 enum EModelType {kONNX, kKeras, kPt, kROOT, kNotDef}; // type of model
62 EModelType type = kNotDef;
63
64 auto pos1 = path.rfind("/");
65 auto pos2 = path.find(".onnx");
66 if (pos2 != std::string::npos) {
67 type = kONNX;
68 } else {
69 pos2 = path.find(".h5");
70 if (pos2 != std::string::npos) {
71 type = kKeras;
72 } else {
73 pos2 = path.find(".pt");
74 if (pos2 != std::string::npos) {
75 type = kPt;
76 }
77 else {
78 pos2 = path.find(".root");
79 if (pos2 != std::string::npos) {
80 type = kROOT;
81 }
82 }
83 }
84 }
85 if (type == kNotDef) {
86 throw std::runtime_error("Input file is not an ONNX or Keras or PyTorch file");
87 }
88 if (pos1 == std::string::npos)
89 pos1 = 0;
90 else
91 pos1 += 1;
92 std::string modelName = path.substr(pos1,pos2-pos1);
93 std::string fileType = path.substr(pos2+1, path.length()-pos2-1);
94 if (verbose) std::cout << "Parsing SOFIE model " << modelName << " of type " << fileType << std::endl;
95
96 // create code for parsing model and generate C++ code for inference
97 // make it in a separate scope to avoid polluting global interpreter space
98 std::string parserCode;
99 if (type == kONNX) {
100 // check first if we can load the SOFIE parser library
101 if (gSystem->Load("libROOTTMVASofieParser") < 0) {
102 throw std::runtime_error("RSofieReader: cannot use SOFIE with ONNX since libROOTTMVASofieParser is missing");
103 }
104 gInterpreter->Declare("#include \"TMVA/RModelParser_ONNX.hxx\"");
105 parserCode += "{\nTMVA::Experimental::SOFIE::RModelParser_ONNX parser ; \n";
106 if (verbose == 2)
107 parserCode += "TMVA::Experimental::SOFIE::RModel model = parser.Parse(\"" + path + "\",true); \n";
108 else
109 parserCode += "TMVA::Experimental::SOFIE::RModel model = parser.Parse(\"" + path + "\"); \n";
110 }
111 else if (type == kKeras) {
112 // use Keras direct parser
113 if (gSystem->Load("libPyMVA") < 0) {
114 throw std::runtime_error("RSofieReader: cannot use SOFIE with Keras since libPyMVA is missing");
115 }
116 // assume batch size is first entry in first input !
117 std::string batch_size = "-1";
118 if (!inputShapes.empty() && ! inputShapes[0].empty())
119 batch_size = std::to_string(inputShapes[0][0]);
120 parserCode += "{\nTMVA::Experimental::SOFIE::RModel model = TMVA::Experimental::SOFIE::PyKeras::Parse(\"" + path +
121 "\"," + batch_size + "); \n";
122 }
123 else if (type == kPt) {
124 // use PyTorch direct parser
125 if (gSystem->Load("libPyMVA") < 0) {
126 throw std::runtime_error("RSofieReader: cannot use SOFIE with PyTorch since libPyMVA is missing");
127 }
128 if (inputShapes.size() == 0) {
129 throw std::runtime_error("RSofieReader: cannot use SOFIE with PyTorch since the input tensor shape is missing and is needed by the PyTorch parser");
130 }
131 std::string inputShapesStr = "{";
132 for (unsigned int i = 0; i < inputShapes.size(); i++) {
133 inputShapesStr += "{ ";
134 for (unsigned int j = 0; j < inputShapes[i].size(); j++) {
135 inputShapesStr += ROOT::Math::Util::ToString(inputShapes[i][j]);
136 if (j < inputShapes[i].size()-1) inputShapesStr += ", ";
137 }
138 inputShapesStr += "}";
139 if (i < inputShapes.size()-1) inputShapesStr += ", ";
140 }
141 inputShapesStr += "}";
142 parserCode += "{\nTMVA::Experimental::SOFIE::RModel model = TMVA::Experimental::SOFIE::PyTorch::Parse(\"" + path + "\", "
143 + inputShapesStr + "); \n";
144 }
145 else if (type == kROOT) {
146 // use parser from ROOT
147 parserCode += "{\nauto fileRead = TFile::Open(\"" + path + "\",\"READ\");\n";
148 parserCode += "TMVA::Experimental::SOFIE::RModel * modelPtr;\n";
149 parserCode += "auto keyList = fileRead->GetListOfKeys(); TString name;\n";
150 parserCode += "for (const auto&& k : *keyList) { \n";
151 parserCode += " TString cname = ((TKey*)k)->GetClassName(); if (cname==\"TMVA::Experimental::SOFIE::RModel\") name = k->GetName(); }\n";
152 parserCode += "fileRead->GetObject(name,modelPtr); fileRead->Close(); delete fileRead;\n";
153 parserCode += "TMVA::Experimental::SOFIE::RModel & model = *modelPtr;\n";
154 }
155
156 // add custom operators if needed
157 if (fCustomOperators.size() > 0) {
158
159 for (auto & op : fCustomOperators) {
160 parserCode += "{ auto p = new TMVA::Experimental::SOFIE::ROperator_Custom<float>(\""
161 + op.fOpName + "\"," + op.fInputNames + "," + op.fOutputNames + "," + op.fOutputShapes + ",\"" + op.fFileName + "\");\n";
162 parserCode += "std::unique_ptr<TMVA::Experimental::SOFIE::ROperator> op(p);\n";
163 parserCode += "model.AddOperator(std::move(op));\n}\n";
164 }
165 }
166
167 int batchSize = 1;
168 if (inputShapes.size() > 0 && inputShapes[0].size() > 0) {
169 batchSize = inputShapes[0][0];
170 if (batchSize < 1) batchSize = 1;
171 }
172 if (verbose) std::cout << "generating the code with batch size = " << batchSize << " ...\n";
173
174 parserCode += "model.Generate(TMVA::Experimental::SOFIE::Options::kDefault,"
175 + ROOT::Math::Util::ToString(batchSize) + ", 0, " + std::to_string(verbose) + "); \n";
176
177 if (verbose) {
178 parserCode += "model.PrintRequiredInputTensors();\n";
179 parserCode += "model.PrintIntermediateTensors();\n";
180 parserCode += "model.PrintOutputTensors();\n";
181 }
182
183 // add custom operators if needed
184#if 0
185 if (fCustomOperators.size() > 0) {
186 if (verbose) {
187 parserCode += "model.PrintRequiredInputTensors();\n";
188 parserCode += "model.PrintIntermediateTensors();\n";
189 parserCode += "model.PrintOutputTensors();\n";
190 }
191 for (auto & op : fCustomOperators) {
192 parserCode += "{ auto p = new TMVA::Experimental::SOFIE::ROperator_Custom<float>(\""
193 + op.fOpName + "\"," + op.fInputNames + "," + op.fOutputNames + "," + op.fOutputShapes + ",\"" + op.fFileName + "\");\n";
194 parserCode += "std::unique_ptr<TMVA::Experimental::SOFIE::ROperator> op(p);\n";
195 parserCode += "model.AddOperator(std::move(op));\n}\n";
196 }
197 parserCode += "model.Generate(TMVA::Experimental::SOFIE::Options::kDefault,"
198 + ROOT::Math::Util::ToString(batchSize) + "); \n";
199 }
200#endif
201 if (verbose > 1)
202 parserCode += "model.PrintGenerated(); \n";
203 parserCode += "model.OutputGenerated();\n";
204
205 parserCode += "int nInputs = model.GetInputTensorNames().size();\n";
206
207 // need information on number of inputs (assume output is 1)
208
209 //end of parsing code, close the scope and return 1 to indicate a success
210 parserCode += "return nInputs;\n}\n";
211
212 if (verbose) std::cout << "//ParserCode being executed:\n" << parserCode << std::endl;
213
214 auto iret = gROOT->ProcessLine(parserCode.c_str());
215 if (iret <= 0) {
216 std::string msg = "RSofieReader: error processing the parser code: \n" + parserCode;
217 throw std::runtime_error(msg);
218 }
219 fNInputs = iret;
220 if (fNInputs > 3) {
221 throw std::runtime_error("RSofieReader does not yet support model with > 3 inputs");
222 }
223
224 // compile now the generated code and create Session class
225 std::string modelHeader = modelName + ".hxx";
226 if (verbose) std::cout << "compile generated code from file " <<modelHeader << std::endl;
227 if (gSystem->AccessPathName(modelHeader.c_str())) {
228 std::string msg = "RSofieReader: input header file " + modelHeader + " is not existing";
229 throw std::runtime_error(msg);
230 }
231 if (verbose) std::cout << "Creating Inference function for model " << modelName << std::endl;
232 std::string declCode;
233 declCode += "#pragma cling optimize(2)\n";
234 declCode += "#include \"" + modelHeader + "\"\n";
235 // create global session instance: use UUID to have an unique name
236 std::string sessionClassName = "TMVA_SOFIE_" + modelName + "::Session";
237 TUUID uuid;
238 std::string uidName = uuid.AsString();
239 uidName.erase(std::remove_if(uidName.begin(), uidName.end(),
240 []( char const& c ) -> bool { return !std::isalnum(c); } ), uidName.end());
241
242 std::string sessionName = "session_" + uidName;
243 declCode += sessionClassName + " " + sessionName + ";";
244
245 if (verbose) std::cout << "//global session declaration\n" << declCode << std::endl;
246
247 bool ret = gInterpreter->Declare(declCode.c_str());
248 if (!ret) {
249 std::string msg = "RSofieReader: error compiling inference code and creating session class\n" + declCode;
250 throw std::runtime_error(msg);
251 }
252
253 fSessionPtr = (void *) gInterpreter->Calc(sessionName.c_str());
254
255 // define a function to be called for inference
256 std::stringstream ifuncCode;
257 std::string funcName = "SofieInference_" + uidName;
258 ifuncCode << "std::vector<float> " + funcName + "( void * ptr";
259 for (int i = 0; i < fNInputs; i++)
260 ifuncCode << ", float * data" << i;
261 ifuncCode << ") {\n";
262 ifuncCode << " " << sessionClassName << " * s = " << "(" << sessionClassName << "*) (ptr);\n";
263 ifuncCode << " return s->infer(";
264 for (int i = 0; i < fNInputs; i++) {
265 if (i>0) ifuncCode << ",";
266 ifuncCode << "data" << i;
267 }
268 ifuncCode << ");\n";
269 ifuncCode << "}\n";
270
271 if (verbose) std::cout << "//Inference function code using global session instance\n"
272 << ifuncCode.str() << std::endl;
273
274 ret = gInterpreter->Declare(ifuncCode.str().c_str());
275 if (!ret) {
276 std::string msg = "RSofieReader: error compiling inference function\n" + ifuncCode.str();
277 throw std::runtime_error(msg);
278 }
279 fFuncPtr = (void *) gInterpreter->Calc(funcName.c_str());
280 //fFuncPtr = reinterpret_cast<std::vector<float> (*)(void *, const float *)>(fptr);
281 fInitialized = true;
282 }
283
284 // Add custom operator
285 void AddCustomOperator(const std::string &opName, const std::string &inputNames, const std::string & outputNames,
286 const std::string & outputShapes, const std::string & fileName) {
287 if (fInitialized) std::cout << "WARNING: Model is already loaded and initialised. It must be done after adding the custom operators" << std::endl;
288 fCustomOperators.push_back( {fileName, opName,inputNames, outputNames,outputShapes});
289 }
290
291 // implementations for different outputs
292 std::vector<float> DoCompute(const std::vector<float> & x1) {
293 if (fNInputs != 1) {
294 std::string msg = "Wrong number of inputs - model requires " + std::to_string(fNInputs);
295 throw std::runtime_error(msg);
296 }
297 auto fptr = reinterpret_cast<std::vector<float> (*)(void *, const float *)>(fFuncPtr);
298 return fptr(fSessionPtr, x1.data());
299 }
300 std::vector<float> DoCompute(const std::vector<float> & x1, const std::vector<float> & x2) {
301 if (fNInputs != 2) {
302 std::string msg = "Wrong number of inputs - model requires " + std::to_string(fNInputs);
303 throw std::runtime_error(msg);
304 }
305 auto fptr = reinterpret_cast<std::vector<float> (*)(void *, const float *, const float *)>(fFuncPtr);
306 return fptr(fSessionPtr, x1.data(),x2.data());
307 }
308 std::vector<float> DoCompute(const std::vector<float> & x1, const std::vector<float> & x2, const std::vector<float> & x3) {
309 if (fNInputs != 3) {
310 std::string msg = "Wrong number of inputs - model requires " + std::to_string(fNInputs);
311 throw std::runtime_error(msg);
312 }
313 auto fptr = reinterpret_cast<std::vector<float> (*)(void *, const float *, const float *, const float *)>(fFuncPtr);
314 return fptr(fSessionPtr, x1.data(),x2.data(),x3.data());
315 }
316
317 /// Compute model prediction on vector
318 template<typename... T>
319 std::vector<float> Compute(T... x)
320 {
321 if(!fInitialized) {
322 return std::vector<float>();
323 }
324
325 // Take lock to protect model evaluation
327
328 // Evaluate TMVA model (need to add support for multiple outputs)
329 return DoCompute(x...);
330
331 }
332 std::vector<float> Compute(const std::vector<float> &x) {
333 if(!fInitialized) {
334 return std::vector<float>();
335 }
336
337 // Take lock to protect model evaluation
339
340 // Evaluate TMVA model (need to add support for multiple outputs)
341 return DoCompute(x);
342 }
343 /// Compute model prediction on input RTensor
344 /// The shape of the input tensor should be {nevents, nfeatures}
345 /// and the return shape will be {nevents, noutputs}
346 /// support for now only a single input
348 {
349 if(!fInitialized) {
350 return RTensor<float>({0});
351 }
352 const auto nrows = x.GetShape()[0];
353 const auto rowsize = x.GetStrides()[0];
354 auto fptr = reinterpret_cast<std::vector<float> (*)(void *, const float *)>(fFuncPtr);
355 auto result = fptr(fSessionPtr, x.GetData());
356
357 RTensor<float> y({nrows, result.size()}, MemoryLayout::ColumnMajor);
358 std::copy(result.begin(),result.end(), y.GetData());
359 //const bool layout = x.GetMemoryLayout() == MemoryLayout::ColumnMajor ? false : true;
360 // assume column major layout
361 for (size_t i = 1; i < nrows; i++) {
362 result = fptr(fSessionPtr, x.GetData() + i*rowsize);
363 std::copy(result.begin(),result.end(), y.GetData() + i*result.size());
364 }
365 return y;
366 }
367
368private:
369
370 bool fInitialized = false;
371 int fNInputs = 0;
372 void * fSessionPtr = nullptr;
373 void * fFuncPtr = nullptr;
374
375 // data to insert custom operators
377 std::string fFileName; // code implementing the custom operator
378 std::string fOpName; // operator name
379 std::string fInputNames; // input tensor names (convert as string as {"n1", "n2"})
380 std::string fOutputNames; // output tensor names converted as trind
381 std::string fOutputShapes; // output shapes
382 };
383 std::vector<CustomOperatorData> fCustomOperators;
384
385};
386
387} // namespace Experimental
388} // namespace TMVA
389
390#endif // TMVA_RREADER
#define c(i)
Definition RSha256.hxx:101
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
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
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
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 type
#define gInterpreter
#define gROOT
Definition TROOT.h:406
R__EXTERN TSystem * gSystem
Definition TSystem.h:561
#define R__WRITE_LOCKGUARD(mutex)
TMVA::RSofieReader class for reading external Machine Learning models in ONNX files,...
RSofieReader(const std::string &path, std::vector< std::vector< size_t > > inputShapes={}, int verbose=0)
Create TMVA model from ONNX file print level can be 0 (minimal) 1 with info , 2 with all ONNX parsing...
RTensor< float > Compute(RTensor< float > &x)
Compute model prediction on input RTensor The shape of the input tensor should be {nevents,...
std::vector< float > Compute(const std::vector< float > &x)
std::vector< float > Compute(T... x)
Compute model prediction on vector.
void Load(const std::string &path, std::vector< std::vector< size_t > > inputShapes={}, int verbose=0)
std::vector< float > DoCompute(const std::vector< float > &x1, const std::vector< float > &x2, const std::vector< float > &x3)
std::vector< CustomOperatorData > fCustomOperators
std::vector< float > DoCompute(const std::vector< float > &x1)
void AddCustomOperator(const std::string &opName, const std::string &inputNames, const std::string &outputNames, const std::string &outputShapes, const std::string &fileName)
std::vector< float > DoCompute(const std::vector< float > &x1, const std::vector< float > &x2)
RSofieReader()
Dummy constructor which needs model loading afterwards.
RTensor is a container with contiguous memory and shape information.
Definition RTensor.hxx:162
const Shape_t & GetShape() const
Definition RTensor.hxx:242
virtual int Load(const char *module, const char *entry="", Bool_t system=kFALSE)
Load a shared library.
Definition TSystem.cxx:1857
virtual Bool_t AccessPathName(const char *path, EAccessMode mode=kFileExists)
Returns FALSE if one can access a file using the specified access mode.
Definition TSystem.cxx:1296
This class defines a UUID (Universally Unique IDentifier), also known as GUIDs (Globally Unique IDent...
Definition TUUID.h:42
const char * AsString() const
Return UUID as string. Copy string immediately since it will be reused.
Definition TUUID.cxx:571
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
std::string ToString(const T &val)
Utility function for conversion to strings.
Definition Util.h:50
R__EXTERN TVirtualRWMutex * gCoreMutex
create variable transformations