Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
mathcoreVectorFloatIO.C File Reference

Detailed Description

View in nbviewer Open in SWAN
Macro illustrating I/O with Lorentz Vectors of floats The dictionary for LorentzVector of float is not in the libMathCore, therefore is generated when parsed the file with CLING.

To run this macro you must do

root[0] .L mathcoreVectorFloatIO.C+
root[1] runIt();
#include "TRandom.h"
#include "TStopwatch.h"
#include "TSystem.h"
#include "TFile.h"
#include "TTree.h"
#include "TH1D.h"
#include "TCanvas.h"
#include <iostream>
#include "TLorentzVector.h"
#include "Math/Vector4D.h"
// Now the dictionary contains the vector's with float types
// No need to force dictionary generation
// You need to run ACLIC with old ROOT version
// and uncomment these lines below
// #ifdef __MAKECINT__
// #pragma link C++ class ROOT::Math::PxPyPzE4D<float>+;
// #pragma link C++ class ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<float> >+;
// #pragma link C++ typedef ROOT::Math::XYZTVectorF;
// #endif
using namespace ROOT::Math;
void write(int n) {
TStopwatch timer;
TFile f1("mathcoreVectorIO_F.root","RECREATE");
// create tree
TTree t1("t1","Tree with new Float LorentzVector");
t1.Branch("LV branch","ROOT::Math::XYZTVectorF",&v1);
timer.Start();
for (int i = 0; i < n; ++i) {
double Px = R.Gaus(0,10);
double Py = R.Gaus(0,10);
double Pz = R.Gaus(0,10);
double E = R.Gaus(100,10);
v1->SetCoordinates(Px,Py,Pz,E);
t1.Fill();
}
f1.Write();
timer.Stop();
std::cout << " Time for new Float Vector " << timer.RealTime() << " " << timer.CpuTime() << std::endl;
t1.Print();
}
void read() {
TStopwatch timer;
TFile f1("mathcoreVectorIO_F.root");
// create tree
TTree *t1 = (TTree*)f1.Get("t1");
t1->SetBranchAddress("LV branch",&v1);
timer.Start();
int n = (int) t1->GetEntries();
std::cout << " Tree Entries " << n << std::endl;
double etot=0;
for (int i = 0; i < n; ++i) {
t1->GetEntry(i);
etot += v1->E();
}
timer.Stop();
std::cout << " Time for new Float Vector " << timer.RealTime() << " " << timer.CpuTime() << std::endl;
std::cout << " E average" << n<< " " << etot << " " << etot/double(n) << endl;
}
void runIt() {
#if defined(__CINT__) && !defined(__MAKECINT__)
gSystem->Load("libMathCore");
gSystem->Load("libPhysics");
using namespace ROOT::Math ;
cout << "This tutorial can run only using ACliC, you must run it by doing: " << endl;
cout << "\t .L tutorials/math/mathcoreVectorFloatIO.C+" << endl;
cout << "\t runIt()" << endl;
#endif
int nEvents = 100000;
write(nEvents);
read();
}
void mathcoreVectorFloatIO() {
#if defined(__CINT__) && !defined(__MAKECINT__)
gSystem->Load("libMathCore");
gSystem->Load("libPhysics");
using namespace ROOT::Math ;
cout << "This tutorial can run only using ACliC, you must run it by doing: " << endl;
cout << "\t .L tutorials/math/mathcoreVectorFloatIO.C+" << endl;
cout << "\t runIt()" << endl;
#endif
}
#define R(a, b, c, d, e, f, g, h, i)
Definition RSha256.hxx:110
R__EXTERN TSystem * gSystem
Definition TSystem.h:555
Class describing a generic LorentzVector in the 4D space-time, using the specified coordinate system ...
A ROOT file is composed of a header, followed by consecutive data records (TKey instances) with a wel...
Definition TFile.h:53
virtual Int_t Write(const char *name=nullptr, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
Definition TObject.cxx:880
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
Stopwatch class.
Definition TStopwatch.h:28
Double_t RealTime()
Stop the stopwatch (if it is running) and return the realtime (in seconds) passed between the start a...
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Double_t CpuTime()
Stop the stopwatch (if it is running) and return the cputime (in seconds) passed between the start an...
void Stop()
Stop the stopwatch.
virtual int Load(const char *module, const char *entry="", Bool_t system=kFALSE)
Load a shared library.
Definition TSystem.cxx:1857
A TTree represents a columnar dataset.
Definition TTree.h:79
const Int_t n
Definition legend1.C:16
TF1 * f1
Definition legend1.C:11
auto * t1
Definition textangle.C:20
Author
Lorenzo Moneta

Definition in file mathcoreVectorFloatIO.C.