STL support

ROOT
1

Hi Rooters,
I would need to fill a TTree with stl:vector because I would like don’t fix their size. I created a class which contains the stl vectors.

[code]#ifndef BDSRoot_h
#define BDSRoot_h
#include
#include “TFile.h”
#include “TTree.h”
#include “TBrowser.h”
#include “TH2.h”
#include “TMath.h”
#include “TROOT.h”
#include “TRandom.h”
#include “TCanvas.h”
#include “TObject.h”
#include

using std::vector;
class BDSRoot : public TObject {
public:

Int_t Nsampler;
Int_t n_hits[200];
std::vector <Float_t> X[200];
std::vector <Float_t> Z[200];;
void Init()
    {
        for(Int_t i=0;i<200;i++)
        {
            X[i].clear();
            Z[i].clear();
        }
    }

void GetPos(Int_t ksampler,Int_t nhit,Float_t x,Float_t z)
    {
        std::cout<<ksampler<<"\t"<<
            nhit<<" x= "<<x<<std::endl;
        
        X[ksampler].push_back(x);
        Z[ksampler].push_back(z);
    }


BDSRoot();
~BDSRoot();

ClassDef(BDSRoot,1)

};
#endif
[/code]

And the code which fills the Tree :

[code]#include “BDSRoot.hh”
#include “BDSSamplerSD.hh”
#include

extern G4String outputFilename;

BDSRoot *hevent=0;

BDSOutput::BDSOutput()
{
//time_t tm = time(NULL);

format = _ASCII; // default - write an ascii file
nSamplers = 0;
//of.open(“output.txt”);
//of<<"## BDSIM output created “<<ctime(&tm)<<” ####"<<G4endl;
}

BDSOutput::BDSOutput(int fmt)
{
format = fmt;
nSamplers = 0;
}

BDSOutput::~BDSOutput()
{
if(of != NULL)
of.close();

#ifdef USE_ROOT
if(format==_ROOT)
if(theRootOutputFile->IsOpen())
{
theRootOutputFile->Write();
theRootOutputFile->Close();
delete theRootOutputFile;
}
#endif
}

void BDSOutput::SetFormat(G4int val)
{
format = val;
time_t tm = time(NULL);

if( format == _ASCII)
{
G4String filename = outputFilename+".txt";
G4cout<<“output format ASCII, filename: “<<filename<<G4endl;
of.open(filename);
of<<”### BDSIM output created “<<ctime(&tm)<<” ####”<<G4endl;
of<<"# PT E[GeV] X[mum] Y[mum] Z[m] Xp[rad] Yp[rad] "<<G4endl;
}
if( format == _ROOT)
{
G4cout<<“output format ROOT”<<G4endl;
}
}

// output initialization (ROOT only)
void BDSOutput::Init(G4int FileNum)
{
if(format != _ROOT) return;

#ifdef USE_ROOT
// set up the root file
G4String filename = outputFilename + “_”
+ BDSGlobals->StringFromInt(FileNum) + “.root”;

G4cout<<"Setting up new file: "<<filename<<G4endl;
theRootOutputFile=new TFile(filename,“RECREATE”, “BDS output file”);

// TTree* SamplerTree = new TTree(“AllDetectors”, “Sampler output”);

TTree* SamplerTree = new TTree(“AllDetectors”, “Sampler output”);

gSystem->Load("/cern/BDSIM_last/rootlib/BDSRoot_cc.so");

G4cout<<"---------->>> DEFINITION DE hevent -------->>>"<<G4endl;
if(!hevent) hevent = new BDSRoot();

// BDSRoot *hevent = new BDSRoot();

int bufsize = 10000;
int split = 1;

SamplerTree->Branch(“hevent”,&hevent, bufsize,split);

hevent->Init();

Nsampler=nSamplers;

for(G4int i=0;i<nSamplers;i++)
{
G4String name=SampName[i];
names[name.c_str()]=i;

  G4cout<<"i = "<<i<<" Name --> "<<name<<G4endl;

}
G4cout<<"----> Apres Declaration des Branches …"<<G4endl;

if(
BDSGlobals->GetStoreTrajectory() ||
BDSGlobals->GetStoreMuonTrajectories() ||
BDSGlobals->GetStoreNeutronTrajectories()
)
// create a tree with trajectories
{
//G4cout<<“BDSOutput::storing trajectories set”<<G4endl;
TTree* TrajTree = new TTree(“Trajectories”, “Trajectories”);
}

// build energy loss histogram
G4int nBins = G4int(zMax/m);

EnergyLossHisto = new TH1F(“ElossHisto”, “Energy Loss”,nBins,0.,zMax/m);
EnergyLossNtuple= new TNtuple(“ElossNtuple”, “Energy Loss”,“z:E”);

#endif
}

void BDSOutput::WriteHits(BDSSamplerHitsCollection *hc)
{
if( format == _ASCII) {

G4cout.precision(6);

for (G4int i=0; i<hc->entries(); i++)
  {
of<<(*hc)[i]->GetPDGtype()
  <<" "
  <<(*hc)[i]->GetMom()/GeV
  <<" "
  <<(*hc)[i]->GetX()/micrometer
  <<" "
  <<(*hc)[i]->GetY()/micrometer
  <<" "
  <<(*hc)[i]->GetZ() / m
  <<" "
  <<(*hc)[i]->GetXPrime() / radian
  <<" "
  <<(*hc)[i]->GetYPrime() / radian
  <<G4endl;
  }

of<<G4endl;
//of<<"end of hits collection"<<G4endl;

}

if( format == _ROOT) {
#ifdef USE_ROOT

// ---- Initialize n_hits ------
for(G4int ks=0; ks<Nsampler;ks++)
{
n_hits[ks]=0;
}

  G4String name;
  
  TTree* sTree=(TTree*)gDirectory->Get("AllDetectors");      
  for (G4int i=0; i<hc->entries(); i++)
  {
      
      G4int ksampler = names[(*hc)[i]->GetName()];
      name=(*hc)[i]->GetName();

      G4cout<<ksampler<<" ICI : "<<i<<"\t"<<n_hits[ksampler]<<G4endl;
      n_hits[ksampler]++;
      
      hevent->GetPos(
          ksampler,
          n_hits[ksampler],
          (*hc)[i]->GetX()/ micrometer,
          (*hc)[i]->GetZ() / m);               

      G4cout<<i<<"before --> "<<ksampler<<"\t"<<
          (*hc)[i]->GetX()/ micrometer<<
          " name --> "<<name.c_str()<< G4endl;
      
      Nsampler=nSamplers;
 }
sTree->Fill();


G4cout<<" Nentries :  "<<hc->entries()<<G4endl;
G4cout<<"------------------------------- "<<G4endl;

#endif
}

}

// write a trajectory to a root/ascii file
G4int BDSOutput::WriteTrajectory(TrajectoryVector* TrajVec)
{
G4cout<<“a trajectory stored…”<<G4endl;

G4int tID;
G4TrajectoryPoint* TrajPoint;
G4ThreeVector TrajPos;

#ifdef USE_ROOT
TTree* TrajTree;

G4String name = “Trajectories”;

TrajTree=(TTree*)gDirectory->Get(name);

if(TrajTree == NULL) { G4cerr<<“TrajTree=NULL”<<G4endl; return -1;}
#endif

if(TrajVec)
{
TrajectoryVector::iterator iT;

  for(iT=TrajVec->begin();iT<TrajVec->end();iT++)
{
  G4Trajectory* Traj=(G4Trajectory*)(*iT);
  
  tID=Traj->GetTrackID();

// part = Traj->GetPDGEncoding();

  for(G4int j=0; j<Traj->GetPointEntries(); j++)
    {
      TrajPoint=(G4TrajectoryPoint*)Traj->GetPoint(j);
      TrajPos=TrajPoint->GetPosition();

// x = TrajPos.x() / m;
// y = TrajPos.y() / m;
// z = TrajPos.z() / m;

#ifdef USE_ROOT
if( format == _ROOT)
TrajTree->Fill();
#endif

      //G4cout<<"TrajPos="<<TrajPos<<G4endl;
    }
}
}

return 0;
}

// make energy loss histo

void BDSOutput::WriteEnergyLoss(BDSEnergyCounterHitsCollection* hc)
{

if( format == _ROOT) {
#ifdef USE_ROOT

G4int n_hit = hc->entries();

for (G4int i=0;i<n_hit;i++)
  {
G4double Energy=(*hc)[i]->GetEnergy();
G4double EWeightZ=(*hc)[i]->
  GetEnergyWeightedPosition()/Energy;

EnergyLossHisto->Fill(EWeightZ/m,Energy/GeV);

EnergyLossNtuple->Fill(EWeightZ/m,Energy/GeV);

  }

#endif
}

if( format == _ASCII) {

G4int n_hit = hc->entries();

for (G4int i=0;i<n_hit;i++)
  {
G4double Energy=(*hc)[i]->GetEnergy();
G4double EWeightZ=(*hc)[i]->
  GetEnergyWeightedPosition()/Energy;

of<<EWeightZ/m<<"  "<<Energy/GeV<<G4endl;

  }

}

}

// write some comments to the output file
// only for ASCII output
void BDSOutput::Echo(G4String str)
{
if(format == _ASCII) of<<"#"<<str<<G4endl;
else // default
G4cout<<"#"<<str<<G4endl;
}

G4int BDSOutput::Commit(G4int FileNum)
{
#ifdef USE_ROOT
if(format==_ROOT)
if(theRootOutputFile->IsOpen())
{
theRootOutputFile->Write();
theRootOutputFile->Close();
delete theRootOutputFile;
}
#endif

Init(FileNum);
return 0;
}
[/code]

So it seems that using that method, the Tree is not properly filled :frowning:
But, if I use instead Float_t X[200][10] it works fine …
And how can I also include a char *name[200] in my Tree ?

Thanks !
Cheers
Hayg

2

Hi,

ROOT does not support array of std::vector and does not support array of char*. Instead use (for example) a std::vector<std::vector > and a std::vector.

Cheers,
Philippe.

3

Thanks for your answer, it works fine !
Except for the TString case … :frowning:
When I fill the rootfile, there is no error in the execution, but as soon as I want to plot it I have the following answer :
Warning in TSelectorDraw::ProcessFillObject: Not implemented for vector
What’s wrong there ?

see the part of the code below :

[quote] std::vector Samplername;
void GetSamplerName(std::string name_tmp)
{
Samplername.push_back(name_tmp);
}[/quote]
Thanks again for the Help !

Cheers,
Hayg

4

Hi,

I can not reproduce this problem.
However, I notice that for a vector of TString, I must explicitly call Data in order to get a real plot:mytree->Scan("obj.myvec.Data()");
while for a single stringmytree->Scan("obj.myTString"); the conversion is automatic.

Cheers,
Philippe.

5

Thanks a lot,
now it works also fine :slight_smile:

Cheers,
Hayg

6

I have a problem concerning the produced rootfile, problem I found recently : I fill the tree like :

[code]std::vector<std::vector<Float_t> > X;
std::vector<std::vector<Float_t> > Y;
std::vector <Float_t> X_temp[200];
std::vector <Float_t> Y_temp[200];

// First : fill the std vectors for each sampler
for (G4int i=0; ientries(); i++)
{ G4int ksampler = names[(*hc)[i]->GetName()];
// from event number, we got the sampler name, and finally its number
X_temp[ksampler].push_back( (X_geant4);
Y_temp[ksampler].push_back( (Y_geant4);
}
// second step : fill the std::vector<std::vector<Float_t>> with the temporary std vectors
for(Int_t is=0;is<Nsampler;is++)
{
X.push_back(X_temp);
Y.push_back(Y_temp);
}
[/code]

The problem comes when I try to plot a distribution from one sampler versus an other sampler, when the number of entries is different in each sampler …
tree->Draw(“X[1]:X[2]”);
Error in TObject::ReadValue: Invalid data address: result will be wrong

So I don’t know from where comes the problem … Is there a real problem for root the handle STL ? Because, in my case, I have samplers with different number of events and I would like to make correlations between them : 1 event in a sampler can correspond to many events in an other sampler.
Thanks a lot for your help !
Cheers, Hayg

7

Hi,

I don’t quite have information to know what’s wrong. Can you provide me with the ROOT file you are trying to use?

Thanks,
Philippe.

8

Thanks for the reply,

See in the Attached file, the TTree :
KEY: TTree AllDetectors;1 Sampler output

===> And please try for example :
AllDetectors->Draw(“X[1]:X[5]”)
—> Error in TObject::ReadValue: Invalid data address: result will be wrong

Cheers,
Hayg
TestSTL.root (809 KB)

9

Hi,

I can indeed reproduce the problem (due to the fact that on some events you have more element in X[1] than in X[5] (or vice et versa) and that this is not handled right.

I will attempt to fix this problem has soon as possible.

Cheers,
Philippe.

10

Hi,

Beside the error message (that I am trying to fix), I notice that you express you intend:[quote]Because, in my case, I have samplers with different number of events and I would like to make correlations between them : 1 event in a sampler can correspond to many events in an other sampler.[/quote]and I am not sure that tree->Draw(“X[1]:X[2]”); when it will work would do that you want. When it works, it will plot the following pairs:

for each i less than min(size(X[1]),size(X[2])) plot( X[1][i] vs X[2][i] );in particular you still have a one to one correspondence (my apologies for not noticing this discrepancy sooner :frowning:).

Cheers,
Philippe.

11

Hi,

The original problem (Error in TObject::ReadValue: Invalid data address: result will be wrong, due to the fact that on some events you have more element in X[1] than in X[5] (or vice et versa)) is fixed in the SVN trunk.

Cheers,
Philippe.

PS. It may still not be doing what you are looking for though :slight_smile:

12

Thanks a lot for providing a new root version which seems to work fine. however, it seems that there is still some problems.
When I plot:

it seems to work fine (Z between 0 and 40) but as soon as I require some conditions on Y[0] like :

the result is quite strange (Z~0 … :unamused: ).
Please use the root file I attached previously to check if the problem is not coming from “my” root version.

Thanks !

13

Use:

For now, setting the index in the condition (by using Y[0]) sets its for the whole expression. So it assumes you wanted to plot Y[0][0]:Z[0].

Cheers,
Philippe.