Hi Helmut,
You are right. Thanks for reporting this problem.
The fix is now in CVS.
Rene Brun
On Fri, 13 Apr 2001, Helmut Schmuecker wrote:
>
> Hi Rene,
>
> Well, the bug is in the first line of:
>
> Bool_t TBits::TestBitNumber(UInt_t bitnumber) const
> {
> // Return the current value of the bit
>
> if (bitnumber>fNbits) return kFALSE;
> UInt_t loc = bitnumber/8;
> UChar_t value = fAllBits[loc];
> UChar_t bit = bitnumber%8;
> Bool_t result = (value & (1<<bit)) != 0;
> return result;
> // short: return 0 != (fAllBits[bitnumber/8] & (1<< (bitnumber%8)));
> }
>
> The ">=" operator must be used there. Otherwise you get an undefined result when
> you ask for bit 8 if it is not yet set (fNbits is set to 8 by the default
> constructor). Sometimes you get a kTRUE, sometimes a kFALSE.
>
> *******************************************
> * *
> * W E L C O M E to R O O T *
> * *
> * Version 3.00/06 10 April 2001 *
> * *
> * You are welcome to visit our Web site *
> * http://root.cern.ch *
> * *
> *******************************************
>
> Compiled with thread support.
>
> CINT/ROOT C/C++ Interpreter version 5.14.79, Feb 24 2001
> Type ? for help. Commands must be C++ statements.
> Enclose multiple statements between { }.
>
> don't panic
>
> root [0] TBits a
> root [1] a.TestBitNumber(8)
> (const unsigned char)0
> root [2] TH1F h1("h1","",1000,1,2)
> root [3] TBits c
> root [4] c.TestBitNumber(8)
> (const unsigned char)1
> root [5]
>
>
> This happens e.g. in
>
> void TFormula::Analyze(const char *schain, Int_t &err)
> ...
> if (k <kMAXFOUND && !fAlreadyFound.TestBitNumber(k)) {
> fAlreadyFound.SetBitNumber(k);
> fNval++;
> }
> ...
>
> Here fNbits of "fAlreadyFound" is set to 8 by the TBits default
> constructor. Thus the result for k==8 is not defined and it may happen that
> fNval is not incremented. fNVal is used in a loop in
>
> Double_t TTreeFormula::EvalInstance(Int_t instance) const
>
> where the necessary leafs for the evaluation of the formula are read from the
> tree. If sombody uses more than 8 leaf names in his selection (TTree::Draw), it
> may happen that the last cut of his selection is applied on some memory trash.
> ... Ouch.
>
> cheers,
> Helmut
>
> long live TCut and TTree::Draw !
> (I have heard that a lot people lost confidence in these fancy methods)
>
>
>
> > Hi Helmut,
> >
> > I cannot reproduce your problem.
> >
> > Rene Brun
> >
> > Helmut Schmuecker wrote:
> > >
> > > Hi rooters,
> > >
> > > After changing from version 2.25-03 to version 3.00-06, one of my macros,
> > > which creates TEventlists by using the TTree::Draw(">>eventlist","some
> > > selection string") method, does not work anymore. With 3.00-06 nothing is
> > > written into the Eventlists if the selection string becomes too long.
> > >
> > > I have tried to reproduce the effect with the TTree T in the
> > > $ROOTSYS/test/Event.root file and the following macro. It is not the same
> > > effect, but what happens is also pretty weird:
> > >
> > > void testtree(){
> > >
> > > gSystem->Load("libEvent.so");
> > > TFile f("Event.root","read");
> > > TTree* tree = (TTree*)f.Get("T");
> > >
> > > TCut cut1 = "abs(fTracks.fBx-0.1)<0.2 && fTracks.fPx>1.0 &&
> > > fTracks.fPz-fTracks.fPy>0.0 && !fTracks.fCharge && !fTracks.fValid &&
> > > fTracks.fPz>1.0 && fTracks.fRandom>800";
> > > TCut cut2 = "fTracks.fNpoint>64";
> > >
> > > TCut magic = cut1 && cut2;
> > >
> > > cerr << "\n\n\n"
> > > << cut1 << "\n\n"
> > > << cut2 << "\n\n"
> > > << magic << endl;
> > >
> > > TCanvas* can = new TCanvas("can","",900,300);
> > > can->Divide(3,1);
> > >
> > > can->cd(1);
> > > T->Draw("fTracks.fNpoint", cut1);
> > > can->cd(2);
> > > T->Draw("fTracks.fNpoint", cut2);
> > > can->cd(3);
> > > T->Draw("fTracks.fNpoint", magic);
> > >
> > > }
> > >
> > > It seems that in "T->Draw("fTracks.fNpoint", magic)" the last cut
> > > (cut2 = "fTracks.fNpoint>64") is simply ignored. In the third histogram
> > > I see also entries with fNpoint<=64. If one of the cuts in
> > > cut1 is removed, e.g. "fTracks.fRandom>800", everything looks fine again.
> > >
> > > With 2.25-03 these problems don't occur.
> > >
> > > Any Idea?
> > >
> > > cheers,
> > > Helmut
> --
> _______________________________________________________________________
>
> Helmut Schmuecker helmut@slac.stanford.edu
> Stanford Linear Accelerator Center helmut@ep1.ruhr-uni-bochum.de
> P.O. Box 4349, MS 41, B 280, R 182 Phone +1-650-926-8653
> Stanford, CA 94309, USA Fax +1-650-926-8522
> _______________________________________________________________________
>
>
>
This archive was generated by hypermail 2b29 : Tue Jan 01 2002 - 17:50:42 MET