Hi,
OK, so I am trying to use ROOT (I thought it would be easier than
writing a Python wrapper for ROOT trees ... turns out I was right,
keep reading ...).
I wrote a macro which ran the very first time (!) which is supposed to
generate a tree full of simulated J/psi electroproduction events.
Each event will contain four TLorentzVector objects. This apparently
happened as far as I can tell.
Pain #1 - reattaching the Tree in the ROOT file seems to be much more
complex than hi/fil 1 charm.root was
Pain #2 - I don't see an easy way to go thru the entire Tree and plot
say the magnitudes of the four-momenta (I thought the tree was
supposed to keep all the object code intact which is why we were using
C++ ... oh dear).
TFile f("charm.root");
TTree *mytree = (TTree*)f.Get("T"); // what a contorted command
mytree.Draw("ppsi.fY"); // yields some error
Error in <TCanvas::Range>: illegal world coordinates range: x1=13.684271, y1=-0.131250, x2=-1.520475, y2=1.181250
Error in <TCanvas::RangeAxis>: illegal axis coordinates range: xmin=12.163796, ymin=0.000000, xmax=-0.000000, ymax=1.050000
mytree.Draw("ppsi.fZ")
Error in <TCanvas::Range>: illegal world coordinates range: x1=13.684271, y1=-0.393750, x2=-1.520475, y2=3.543750
Error in <TCanvas::RangeAxis>: illegal axis coordinates range: xmin=12.163796, ymin=0.000000, xmax=-0.000000, ymax=3.150000
mytree.Draw("ppsi.Mag()")
*ERROR 26 :
Unknown name : "ppsi.Mag()"
Error in <TCanvas::Range>: illegal world coordinates range: x1=13.684271, y1=-0.393750, x2=-1.520475, y2=3.543750
Error in <TCanvas::RangeAxis>: illegal axis coordinates range: xmin=12.163796, ymin=0.000000, xmax=-0.000000, ymax=3.150000
ppsi.Mag() is the biggest disappointment. What is the use of having
all these objects in the tree if I can't use their methods?
Pain #3 - wading thru the documentation gives me no clue of how to do
what I want to do, except for possibly the contorted, non-simple
method shown under "The General Way" on URL
http://root.cern.ch/root/HowtoReadTree.html
Say it ain't so!!! Isn't OO supposed to make my life easier?
JT
ps In the very possible case that I have already blown it in my macro,
here it is.
// #-> charm.py -- simulate j/psi electroproduction at JLab
// c++ version i hope
// from Numeric import * [ left over from python version ]
// from vec4 import *
// import RNG
// #-> simulation driver parameters
main(int argc, char **argv)
{
float mj = 3096.88; // # J/psi mass in MeV
float me = 0.510999; // # electron mass
float mp = 938.272;
float beamen = 11.0 * 1000.0; // # GeV * 1000 MeV/GeV
float ll_the = (3.1415926 / 180) * 1.0; // last number is e- angle LL in deg
float ul_the = (3.1415926 / 180) * 20.0; // last number is e- angle UL in deg
// above is stored in radians since we will need to take cosines
// below is not
float ll_phe = -93.0; // LL on e- azimuthal angle (degrees now)
float ul_phe = -87.0;
float ll_pe = 2275; // LL on scat electron momentum (MeV/c)
float ul_pe = 2780;
int howmany = 10000; // # desired number of events
// #-> set up random generators for electron. strategy: create sample
// #-> arrays at start of program, just take values from them during
// #-> execution
TF1 *thegen = new TF1("thegen","sin(x)",ll_the,ul_the);
TF1 *phegen = new TF1("phegen","1.0", ll_phe,ul_phe);
TF1 *pegen = new TF1("pegen", "1.0", ll_pe, ul_pe );
// #-> set up four-momenta for beam and target, which should not change
TLorentzVector beam; beam.SetVectM(TVector3(0,0,beamen),me);
TLorentzVector targ(0.0,0.0,0.0,mp);
TFile hfile("charm.root","RECREATE","Charm");
TTree *tree = new TTree("T","Charm Electroproduction");
int nev = 0;
TLorentzVector scat(1,1,1,1);
TLorentzVector X;
tree->Branch("beam",&beam,"fX/D:fY/D:fZ/D:fE/D");
tree->Branch("targ",&targ,"fX/D:fY/D:fZ/D:fE/D");
tree->Branch("scat",&scat,"fX/D:fY/D:fZ/D:fE/D");
tree->Branch("ppsi",&X,"fX/D:fY/D:fZ/D:fE/D");
for ( int nev = 1; nev <= howmany; nev++) {
float thtmp = thegen.GetRandom();
float phtmp = phegen.GetRandom();
float ptmp = pegen.GetRandom();
scat.SetVectM(TVector3(ptmp*sin(thtmp)*cos(phtmp),
ptmp*sin(thtmp)*sin(phtmp),
ptmp*cos(thtmp)),
me);
X = beam + targ - scat;
tree->Fill();
}
hfile.Write();
hfile.Close();
return 0;
}
This archive was generated by hypermail 2b29 : Tue Jan 02 2001 - 11:50:17 MET