Re: Writing arrays of TH1F to a file

From: Fons Rademakers (Fons.Rademakers@cern.ch)
Date: Thu Jan 06 2000 - 18:57:02 MET


Hi Norbert,

  actually there is no bug in the rootcint generated code
and for me it works fine. I use:

//------------------------------ MT.h
#ifndef __MT_HH__
#define __MT_HH__
#include <TObject.h>
#include <TH1.h>   // need include, forward decl of TH1F is not enough

class MT: public TObject {

private:
  Float_t FoilZ; // Z Position of the Foil
  TH1F* hSpinAngleCut[256][16]; //!
  TH1F* hSpinAnglePhi;
  TH1F* hSpinAngleDiv[256][16]; //

public:
  MT(){;}
  MT(Text_t* Name, Text_t* Title);
 ~MT(){;}
  void Print(Option_t *);

  ClassDef(MT,1) // Muon Spin Rotation Analysis
};

//--------------------------------- MT.cxx
#include <TString.h>
#include "MT.h"

ClassImp(MT)

MT::MT(Text_t *name, Text_t *title)
{
   for (int i = 0; i < 256; i++)
      for (int j = 0; j < 16; j++)
         hSpinAngleDiv[i][j] = new TH1F(Form("%s%d",name, i*100+j), title,
                                        100,0.,100.);
}

void MT::Print(Option_t *)
{
   for (int i = 0; i < 256; i++)
      for (int j = 0; j < 16; j++)
         hSpinAngleDiv[i][j]->Print();
}

#---------- Make
rootcint -f dict.cxx -c MT.h
g++ -fPIC -g `root-config --cflags` -c MT.cxx
g++ -fPIC -g `root-config --cflags` -c dict.cxx
g++ -g -shared -o libmt.so MT.o dict.o

//---------- write.C
{
   gSystem.Load("libmt");
   MT a("aap","noot");
   TFile *f = new TFile("aa.root","recreate");
   a.Write();
   delete f;
}

//---------- read.C
{
   gSystem.Load("libmt");
   TFile *f = new TFile("aa.root");
   MT.Dump();
   MT.Print();
}


Try this and if it still does not work let me know.

Cheers, Fons.



Rene Brun wrote:
> 
> Hi Norbert,
> The rootcint code generating the Streamer function in case
> of a 2-d array of pointers is not correct.
> We will fix this in the next release.
> Meanwhile implement your owm Streamer.
> 
> Rene Brun
> 
> Norbert Danneberg wrote:
> >
> > HI ,
> >
> > I am writing an Object (MTMuSRAnalysis) which contains a 2-dim  array of
> > TH1F to a file.
> > When I read this object back from disk I get the error :
> >
> > Error in <TBuffer::ReadObject>: got object of wrong class
> > Error in <TBuffer::ReadObject>: got object of wrong class
> >
> > This problem does not exist in the case of one dimensional arrays of
> > histograms.
> >
> > Do I have to write a custum streamer function ?
> >
> > Cheers,
> >
> > Norbert
> >
> > class MTMuSRAnalysis: public MTAnalysis
> > {
> > private:
> >   Float_t FoilZ; // Z Position of the Foil
> >   TH1F* hSpinAngleCut[256][16]; //!
> >   TH1F* hSpinAnglePhi;
> >   TH1F* hSpinAngleDiv[256][16]; //
> >
> > public:
> >   MTMuSRAnalysis(){;}
> >   MTMuSRAnalysis(Text_t* Name, Text_t* Title, TTree* tree);
> >  ~MTMuSRAnalysis(){;}
> >
> >   ClassDef(MTMuSRAnalysis,1) // Muon Spin Rotation Analysis
> > };
> >
> > --
> >
> > Norbert Danneberg
> >
> > ETH Zurich - Institute for Particle Physics
> >   Laboratory for Nuclear Physics          Phone.: +41-1-633-2034
> >   Hoenggerberg                            Fax.:   +41-1-633-1067
> >   CH-8093 Zurich
> >
> > ETH Zurich - Institute for Particle Physics
> >   Paul Scherrer Institute                 Phone.: +41-56-310-3284
> >   CH-5232 Villigen PSI                    Fax.:   +41-56-310-4362
> >
> > email: Norbert.Danneberg@psi.ch

-- 
Org:    CERN, European Laboratory for Particle Physics.
Mail:   1211 Geneve 23, Switzerland
E-Mail: Fons.Rademakers@cern.ch              Phone: +41 22 7679248
WWW:    http://root.cern.ch/~rdm/            Fax:   +41 22 7677910



This archive was generated by hypermail 2b29 : Tue Jan 02 2001 - 11:50:16 MET