Hello,
I have an odd problem. Namely when I try to save my data to a root file,
then for some reason I only get the first entry for every event written to
the ROOT file.
In a bit more details:
1. I used a sample ORCA code called ExTutTree to reconstruct ECAL towers.
I modified the code to also reconstruct the muons.
2. the code works 100% in it's original form and it also does the same job
in the modified version (only additions, no removal of code). So there are
approximately 1100 records written to the ROOT tree for every event
reconstructed.
3. when I added the muon reconstruction then I also added the necessary
code to save a few more branches to the ROOT tree. The code works just
fine in theory as it is saving 1 entry per event and the data it saves
corresponds to that what I'm printing on the screen during reconstruction.
4. the problems come when the system encounters more than one muon per
event and tries to save the second one as well. The code that should write
the data to ROOT file is executed (I added some debug text to it and it
gets printed), but the data never makes it to the ROOT file.
Now to allow you to diagnose the stuff I will try to copy some relevant
parts of the code here. The whole code can be looked at
/afs/cern.ch/user/m/mario/public/Tutorial (it's a copy of the working area
that I have in ORCA so it's not the whole code, but the only relevant part
of it).
Firstly there is a separate class for actions with root trees. It's called
ExSimpleCaloTree (with header in interface folder). The actual analysis is
done in ExSimpleCaloAnalysis and the code in Analysis is executed for
every event.
For those of you who can't access the afs path I will try to copy some of
the code:
from interfaces/ExSimpleCaloTree.h
in public:
void addTower(float e, float phi, float eta);
void addMuons(float px, float py, float pz, float p, int n);
void store();
in private:
static const int NTWMAX = 5000;
static const int NMUMAX = 5000;
int myNumberOfTowers;
float myETw[NTWMAX], myPhiTw[NTWMAX], myEtaTw[NTWMAX];
// MUONS
int myMuonNr, myN[NMUMAX];
float myPx[NMUMAX], myPy[NMUMAX], myPz[NMUMAX], myP[NMUMAX];
from ExSimpleCaloTree.cc
void ExSimpleCaloTree::store()
{
myTree->Fill();
myNumberOfTowers=0;
myMuonNr=0;
}
(in constructor)
myFile = new TFile(filename,"RECREATE");
myTree = new TTree("T1","ORCA tutorial tree");
myNumberOfTowers=0;
myTree->Branch("Ntowers", &myNumberOfTowers, "Ntowers/I");
myTree->Branch("Etower", &myETw, "Etower[Ntowers]/F");
myTree->Branch("Phitower", &myPhiTw, "Phitower[Ntowers]/F");
myTree->Branch("Etatower", &myEtaTw, "Etatower[Ntowers]/F");
myMuonNr=0;
myTree->Branch("MuonNr", &myMuonNr, "MuonNr/I");
myTree->Branch("Px", &myPx, "Px/F");
myTree->Branch("Py", &myPy, "Py/F");
myTree->Branch("Pz", &myPz, "Pz/F");
myTree->Branch("P", &myP, "P/F");
myTree->Branch("N", &myN, "N/I");
(in destructor)
myFile->cd();
myTree->Write();
myFile->Close();
delete myFile;
(actual writing functions)
void ExSimpleCaloTree::addTower(float e, float phi, float eta)
{
if (myNumberOfTowers < NTWMAX) {
myETw[myNumberOfTowers]=e;
myPhiTw[myNumberOfTowers]= (phi < 0.) ? (2*M_PI+phi) : phi;
myEtaTw[myNumberOfTowers]=eta;
myNumberOfTowers++;
}
}
void ExSimpleCaloTree::addMuons(float px, float py, float pz, float p, int
n)
{
cout << "myMuonNr = " << myMuonNr << " and NMUMAX " << NMUMAX <<
endl;
if (myMuonNr < NMUMAX) {
cout << "WE ARE SAVING TO TREE" << endl;
myPx[myMuonNr]=px;
myPy[myMuonNr]=py;
myPz[myMuonNr]=pz;
myP[myMuonNr]=p;
myN[myMuonNr]=n;
myMuonNr++;
cout << "So n is " << myN[myMuonNr-1] << " in Tree..." <<
endl;
}
}
so I should get "WE ARE SAVING TO TREE" every time a muon is saved to the
Tree object. Well I do get that text on console for every muon that is
reconstructed.
Here is the relevant part from ExSimpleCaloAnalysis.cc
(the constructor)
ExSimpleCaloAnalysis::ExSimpleCaloAnalysis() :
theInitialised(false),
theCaloTowers(0), Muons(0) {
cout << "Constructing SimpleCalo Observer" << endl;
Observer< G3EventProxy *>::init();
Observer< G3SetUp *>::init();
myTree = new ExSimpleCaloTree("muons.root"); // instantiate the user
tree
};
// this is the destructor
ExSimpleCaloAnalysis::~ExSimpleCaloAnalysis() {
cout << "Deleting SimpleCalo Observer" << endl;
delete myTree;
delete theCaloTowers;
delete Muons;
};
(here we actually loop over all Calo towers and save them to tree)
if ( theCaloTowers==0 )
theCaloTowers = new
RecCollection<EcalPlusHcalTower>(RecQuery("TowerBuilder","1.0"));
cout << "Number of calo towers = " << theCaloTowers->size() << endl;
for (CaloTowerCollection::const_iterator it=theCaloTowers->begin();
it!=theCaloTowers->end(); it++) {
theEcalo += (*it)->energy(); // sum up the total energy
theEecaltotal+=(*it)->energyEcalTower(); // sum up the Ecal energy
theEhcaltotal+=(*it)->energyHcalTower(); // sum up the Hcal energy
theTowerPosition = (*it)->position();
myTree->addTower((*it)->energy(),theTowerPosition.phi(),
theTowerPosition.pseudoRapidity());
}
(and here we reconstruct and save the muons)
if (Muons == 0)
Muons = new RecCollection<RecMuon>(
RecQuery("StandAloneMuonReconstructor"));
GlobalVector p;
int n=0;
cout << "Collection has " << Muons->size() << " muons" << endl;
for(MuonCollection::const_iterator it=Muons->begin();
it!=Muons->end(); it++) {
FreeTrajectoryState* traj = ( (*it)->vertexState() == 0 ) ?
(*it)->innermostState().freeState() :
(*it)->vertexState();
p = traj->momentum();
n++;
myTree->addMuons(p.x(),p.y(),p.z(),p.mag(),n);
cout << "P magnitude: " << p.mag() << endl;
}
myTree->store();
so as you see there shouldn't be any difference in the calorimetry and in
the muon parts. But the calorimeter saves all the data and the muons get
saved only the first one, no others. I see for every muon that the
addMuons is called as the WE ARE SAVING TO TREE text appears on console,
but for a 200 event run I will get 200 muons although I actually
reconstructed 392.
As I don't seem to figure out how the code distinguishes between the first
muon and all the rest and especially that it doesn't behave like this with
the calo towers is more than my head can currently take. Any hints or
help would be appreciated.
for the record: I'm using ORCA_8_2_0 that uses ROOT 3.10 and for analysis
I use the latest ROOT 4.00/08 but I have also tried 3.10 from the ORCA
tree, no difference.
Thank you in advance,
--
Mario Kadastik
National Institute of Chemical Physics and Biophysics
CMS experiment
Mario.Kadastik@cern.ch
--
"Physics is like sex. Sure, it may give some practical results, but that's
not why we do it."
- R. Feynman
This archive was generated by hypermail 2b29 : Sun Jan 02 2005 - 05:50:08 MET