mathcoreVectorIO.C

Go to the documentation of this file.
00001 //
00002 // Example of  I/O of a mathcore Lorentz Vectors in a Tree and comparison with a TLorentzVector
00003 // A ROOT tree is written and read in both using either a XYZTVector or /// a TLorentzVector. 
00004 // 
00005 //  To execute the macro type in: 
00006 //
00007 // root[0]: .x  mathcoreVectorIO.C
00008 //Author: Lorenzo Moneta
00009 
00010 
00011 #include "TRandom2.h"
00012 #include "TStopwatch.h"
00013 #include "TSystem.h"
00014 #include "TFile.h"
00015 #include "TTree.h"
00016 #include "TH1D.h"
00017 #include "TCanvas.h"
00018 
00019 #include <iostream>
00020 
00021 #include "TLorentzVector.h"
00022 
00023 #include "Math/Vector4D.h"
00024 
00025 using namespace ROOT::Math;
00026 
00027  
00028 
00029 
00030 void write(int n) { 
00031 
00032 
00033   TRandom2 R; 
00034   TStopwatch timer;
00035 
00036 
00037   R.SetSeed(1);
00038   timer.Start();
00039   double s = 0; 
00040   for (int i = 0; i < n; ++i) { 
00041         s  += R.Gaus(0,10);
00042         s  += R.Gaus(0,10);
00043         s  += R.Gaus(0,10);
00044         s  += R.Gaus(100,10);
00045   }
00046 
00047   timer.Stop();
00048   std::cout << s/double(n) << std::endl;
00049   std::cout << " Time for Random gen " << timer.RealTime() << "  " << timer.CpuTime() << std::endl; 
00050   
00051 
00052   TFile f1("mathcoreVectorIO_1.root","RECREATE");
00053 
00054   // create tree
00055   TTree t1("t1","Tree with new LorentzVector");
00056 
00057   XYZTVector *v1 = new XYZTVector();
00058   t1.Branch("LV branch","ROOT::Math::XYZTVector",&v1);
00059 
00060   R.SetSeed(1);
00061   timer.Start();
00062   for (int i = 0; i < n; ++i) { 
00063         double Px = R.Gaus(0,10);
00064         double Py = R.Gaus(0,10);
00065         double Pz = R.Gaus(0,10);
00066         double E  = R.Gaus(100,10);
00067         //CylindricalEta4D<double> & c = v1->Coordinates();
00068         //c.SetValues(Px,pY,pZ,E);
00069         v1->SetCoordinates(Px,Py,Pz,E);
00070         t1.Fill(); 
00071   }
00072 
00073   f1.Write();
00074   timer.Stop();
00075   std::cout << " Time for new Vector " << timer.RealTime() << "  " << timer.CpuTime() << std::endl; 
00076 
00077   t1.Print();
00078 
00079   // create tree with old LV 
00080 
00081   TFile f2("mathcoreVectorIO_2.root","RECREATE");
00082   TTree t2("t2","Tree with TLorentzVector");
00083 
00084   TLorentzVector * v2 = new TLorentzVector();
00085   TLorentzVector::Class()->IgnoreTObjectStreamer();
00086   TVector3::Class()->IgnoreTObjectStreamer();
00087 
00088   t2.Branch("TLV branch","TLorentzVector",&v2,16000,2);
00089 
00090   R.SetSeed(1);
00091   timer.Start();
00092   for (int i = 0; i < n; ++i) { 
00093         double Px = R.Gaus(0,10);
00094         double Py = R.Gaus(0,10);
00095         double Pz = R.Gaus(0,10);
00096         double E  = R.Gaus(100,10);
00097         v2->SetPxPyPzE(Px,Py,Pz,E);
00098         t2.Fill(); 
00099   }
00100 
00101   f2.Write();
00102   timer.Stop();
00103   std::cout << " Time for old Vector " << timer.RealTime() << "  " << timer.CpuTime() << endl; 
00104 
00105   t2.Print();
00106 }
00107 
00108 
00109 void read() { 
00110 
00111 
00112 
00113 
00114   TRandom R; 
00115   TStopwatch timer;
00116 
00117 
00118 
00119   TFile f1("mathcoreVectorIO_1.root");
00120 
00121   // create tree
00122   TTree *t1 = (TTree*)f1.Get("t1");
00123 
00124   XYZTVector *v1 = 0;
00125   t1->SetBranchAddress("LV branch",&v1);
00126 
00127   timer.Start();
00128   int n = (int) t1->GetEntries();
00129   std::cout << " Tree Entries " << n << std::endl; 
00130   double etot=0;
00131   for (int i = 0; i < n; ++i) { 
00132     t1->GetEntry(i);
00133     etot += v1->Px();
00134     etot += v1->Py();
00135     etot += v1->Pz();
00136     etot += v1->E();
00137   }
00138 
00139 
00140   timer.Stop();
00141   std::cout << " Time for new Vector " << timer.RealTime() << "  " << timer.CpuTime() << std::endl; 
00142 
00143   std::cout << " TOT average : n = " << n << "\t " << etot/double(n) << endl; 
00144 
00145 
00146   // create tree with old LV 
00147 
00148   TFile f2("mathcoreVectorIO_2.root");
00149   TTree *t2 = (TTree*)f2.Get("t2");
00150 
00151 
00152   TLorentzVector * v2 = 0;
00153   t2->SetBranchAddress("TLV branch",&v2);
00154 
00155   timer.Start();
00156   n = (int) t2->GetEntries();
00157   std::cout << " Tree Entries " << n << std::endl; 
00158   etot = 0;
00159   for (int i = 0; i < n; ++i) { 
00160     t2->GetEntry(i);
00161     etot  += v2->Px();
00162     etot  += v2->Py();
00163     etot  += v2->Pz();
00164     etot  += v2->E();
00165   }
00166 
00167   timer.Stop();
00168   std::cout << " Time for old Vector " << timer.RealTime() << "  " << timer.CpuTime() << endl; 
00169   std::cout << " TOT average:\t" << etot/double(n) << endl; 
00170 }
00171 
00172 
00173 
00174 void mathcoreVectorIO() { 
00175 
00176 
00177 #if defined(__CINT__) && !defined(__MAKECINT__) 
00178 
00179   gSystem->Load("libMathCore");  
00180   gSystem->Load("libPhysics");  
00181   // in CINT need to do that after having loading the library
00182   using namespace ROOT::Math;
00183 
00184   cout << "This tutorial can run only using ACliC, compiling it by doing: " << endl;
00185   cout << "\t  .x tutorials/math/mathcoreVectorCollection.C+" << endl; 
00186   //gROOT->ProcessLine(".x tutorials/math/mathcoreVectorCollection.C+"); 
00187   return;
00188 
00189 #endif
00190 
00191   
00192   int nEvents = 100000;
00193 
00194   write(nEvents);
00195 
00196   read();
00197 }
00198   

Generated on Tue Jul 5 15:44:50 2011 for ROOT_528-00b_version by  doxygen 1.5.1