RooDataHist.cxx

Go to the documentation of this file.
00001 /*****************************************************************************
00002  * Project: RooFit                                                           *
00003  * Package: RooFitCore                                                       *
00004  * @(#)root/roofitcore:$Id: RooDataHist.cxx 37128 2010-11-30 22:24:50Z wouter $
00005  * Authors:                                                                  *
00006  *   WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu       *
00007  *   DK, David Kirkby,    UC Irvine,         dkirkby@uci.edu                 *
00008  *                                                                           *
00009  * Copyright (c) 2000-2005, Regents of the University of California          *
00010  *                          and Stanford University. All rights reserved.    *
00011  *                                                                           *
00012  * Redistribution and use in source and binary forms,                        *
00013  * with or without modification, are permitted according to the terms        *
00014  * listed in LICENSE (http://roofit.sourceforge.net/license.txt)             *
00015  *****************************************************************************/
00016 
00017 //////////////////////////////////////////////////////////////////////////////
00018 // 
00019 // BEGIN_HTML
00020 // RooDataSet is a container class to hold N-dimensional binned data. Each bins central 
00021 // coordinates in N-dimensional space are represented by a RooArgSet of RooRealVar, RooCategory 
00022 // or RooStringVar objects, thus data can be binned in real and/or discrete dimensions
00023 // END_HTML
00024 //
00025 //
00026 
00027 #include "RooFit.h"
00028 #include "Riostream.h"
00029 
00030 #include "TH1.h"
00031 #include "TH1.h"
00032 #include "TDirectory.h"
00033 #include "TMath.h"
00034 #include "RooMsgService.h"
00035 #include "RooDataHist.h"
00036 #include "RooDataHistSliceIter.h"
00037 #include "RooAbsLValue.h"
00038 #include "RooArgList.h"
00039 #include "RooRealVar.h"
00040 #include "RooMath.h"
00041 #include "RooBinning.h"
00042 #include "RooPlot.h"
00043 #include "RooHistError.h"
00044 #include "RooCategory.h"
00045 #include "RooCmdConfig.h"
00046 #include "RooTreeDataStore.h"
00047 #include "TTree.h"
00048 #include "RooTreeData.h"
00049 
00050 using namespace std ;
00051 
00052 ClassImp(RooDataHist) 
00053 ;
00054 
00055 
00056 
00057 //_____________________________________________________________________________
00058 RooDataHist::RooDataHist() : _pbinvCacheMgr(0,10)
00059 {
00060   // Default constructor
00061   _arrSize = 0 ;
00062   _wgt = 0 ;
00063   _errLo = 0 ;
00064   _errHi = 0 ;
00065   _sumw2 = 0 ;
00066   _binv = 0 ;
00067   _pbinv = 0 ;
00068   _curWeight = 0 ;
00069   _curIndex = -1 ;
00070   _realIter = _realVars.createIterator() ;
00071   _binValid = 0 ;
00072 
00073 }
00074 
00075 
00076 
00077 //_____________________________________________________________________________
00078 RooDataHist::RooDataHist(const char *name, const char *title, const RooArgSet& vars, const char* binningName) : 
00079   RooAbsData(name,title,vars), _wgt(0), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00080 {
00081   // Constructor of an empty data hist from a RooArgSet defining the dimensions
00082   // of the data space. The range and number of bins in each dimensions are taken
00083   // from getMin()getMax(),getBins() of each RooAbsArg representing that
00084   // dimension.
00085   //
00086   // For real dimensions, the fit range and number of bins can be set independently
00087   // of the plot range and number of bins, but it is advisable to keep the
00088   // ratio of the plot bin width and the fit bin width an integer value.
00089   // For category dimensions, the fit ranges always comprises all defined states
00090   // and each state is always has its individual bin
00091   //
00092   // To effective achive binning of real dimensions with variable bins sizes,
00093   // construct a RooThresholdCategory of the real dimension to be binned variably.
00094   // Set the thresholds at the desired bin boundaries, and construct the
00095   // data hist as function of the threshold category instead of the real variable.
00096 
00097   // Initialize datastore
00098   _dstore = new RooTreeDataStore(name,title,_vars) ;
00099   
00100   initialize(binningName) ;
00101 
00102   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00103 
00104   appendToDir(this,kTRUE) ;
00105 }
00106 
00107 
00108 
00109 //_____________________________________________________________________________
00110 RooDataHist::RooDataHist(const char *name, const char *title, const RooArgSet& vars, const RooAbsData& data, Double_t wgt) :
00111   RooAbsData(name,title,vars), _wgt(0), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00112 {
00113   // Constructor of a data hist from an existing data collection (binned or unbinned)
00114   // The RooArgSet 'vars' defines the dimensions of the histogram. 
00115   // The range and number of bins in each dimensions are taken
00116   // from getMin()getMax(),getBins() of each RooAbsArg representing that
00117   // dimension.
00118   //
00119   // For real dimensions, the fit range and number of bins can be set independently
00120   // of the plot range and number of bins, but it is advisable to keep the
00121   // ratio of the plot bin width and the fit bin width an integer value.
00122   // For category dimensions, the fit ranges always comprises all defined states
00123   // and each state is always has its individual bin
00124   //
00125   // To effective achive binning of real dimensions with variable bins sizes,
00126   // construct a RooThresholdCategory of the real dimension to be binned variably.
00127   // Set the thresholds at the desired bin boundaries, and construct the
00128   // data hist as function of the threshold category instead of the real variable.
00129   //
00130   // If the constructed data hist has less dimensions that in source data collection,
00131   // all missing dimensions will be projected.
00132 
00133   // Initialize datastore
00134   _dstore = new RooTreeDataStore(name,title,_vars) ;
00135   initialize() ;
00136   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00137 
00138   add(data,(const RooFormulaVar*)0,wgt) ;
00139   appendToDir(this,kTRUE) ;
00140 }
00141 
00142 
00143 
00144 //_____________________________________________________________________________
00145 RooDataHist::RooDataHist(const char *name, const char *title, const RooArgList& vars, RooCategory& indexCat, 
00146                          map<string,TH1*> histMap, Double_t wgt) :
00147   RooAbsData(name,title,RooArgSet(vars,&indexCat)), 
00148   _wgt(0), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00149 {
00150   // Constructor of a data hist from a map of TH1,TH2 or TH3 that are collated into a x+1 dimensional
00151   // RooDataHist where the added dimension is a category that labels the input source as defined
00152   // in the histMap argument. The state names used in histMap must correspond to predefined states
00153   // 'indexCat'
00154   //
00155   // The RooArgList 'vars' defines the dimensions of the histogram. 
00156   // The ranges and number of bins are taken from the input histogram and must be the same in all histograms
00157 
00158   // Initialize datastore
00159   _dstore = new RooTreeDataStore(name,title,_vars) ; 
00160   
00161   importTH1Set(vars, indexCat, histMap, wgt, kTRUE) ;
00162 
00163   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00164 }
00165 
00166 
00167 
00168 //_____________________________________________________________________________
00169 RooDataHist::RooDataHist(const char *name, const char *title, const RooArgList& vars, RooCategory& indexCat, 
00170                          map<string,RooDataHist*> dhistMap, Double_t wgt) :
00171   RooAbsData(name,title,RooArgSet(vars,&indexCat)), 
00172   _wgt(0), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00173 {
00174   // Constructor of a data hist from a map of RooDataHists that are collated into a x+1 dimensional
00175   // RooDataHist where the added dimension is a category that labels the input source as defined
00176   // in the histMap argument. The state names used in histMap must correspond to predefined states
00177   // 'indexCat'
00178   //
00179   // The RooArgList 'vars' defines the dimensions of the histogram. 
00180   // The ranges and number of bins are taken from the input histogram and must be the same in all histograms
00181 
00182   // Initialize datastore
00183   _dstore = new RooTreeDataStore(name,title,_vars) ; 
00184   
00185   importDHistSet(vars, indexCat, dhistMap, wgt) ;
00186 
00187   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00188 }
00189 
00190 
00191 
00192 //_____________________________________________________________________________
00193 RooDataHist::RooDataHist(const char *name, const char *title, const RooArgList& vars, const TH1* hist, Double_t wgt) :
00194   RooAbsData(name,title,vars), _wgt(0), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00195 {
00196   // Constructor of a data hist from an TH1,TH2 or TH3
00197   // The RooArgSet 'vars' defines the dimensions of the histogram. The ranges
00198   // and number of bins are taken from the input histogram, and the corresponding
00199   // values are set accordingly on the arguments in 'vars'
00200 
00201   // Initialize datastore
00202   _dstore = new RooTreeDataStore(name,title,_vars) ; 
00203 
00204   // Check consistency in number of dimensions
00205   if (vars.getSize() != hist->GetDimension()) {
00206     coutE(InputArguments) << "RooDataHist::ctor(" << GetName() << ") ERROR: dimension of input histogram must match "
00207                           << "number of dimension variables" << endl ;
00208     assert(0) ; 
00209   }
00210 
00211   importTH1(vars,*const_cast<TH1*>(hist),wgt, kTRUE) ;
00212 
00213   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00214 }
00215 
00216 
00217 
00218 //_____________________________________________________________________________
00219 RooDataHist::RooDataHist(const char *name, const char *title, const RooArgList& vars, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
00220                          const RooCmdArg& arg4,const RooCmdArg& arg5,const RooCmdArg& arg6,const RooCmdArg& arg7,const RooCmdArg& arg8) :
00221   RooAbsData(name,title,RooArgSet(vars,(RooAbsArg*)RooCmdConfig::decodeObjOnTheFly("RooDataHist::RooDataHist", "IndexCat",0,0,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8))), 
00222   _wgt(0), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00223 {
00224   // Constructor of a binned dataset from a RooArgSet defining the dimensions
00225   // of the data space. The range and number of bins in each dimensions are taken
00226   // from getMin() getMax(),getBins() of each RooAbsArg representing that
00227   // dimension.
00228   //
00229   // This constructor takes the following optional arguments
00230   //
00231   // Import(TH1&, Bool_t impDens) -- Import contents of the given TH1/2/3 into this binned dataset. The 
00232   //                                 ranges and binning of the binned dataset are automatically adjusted to
00233   //                                 match those of the imported histogram. 
00234   //
00235   //                                 Please note: for TH1& with unequal binning _only_,
00236   //                                 you should decide if you want to import the absolute bin content,
00237   //                                 or the bin content expressed as density. The latter is default and will
00238   //                                 result in the same histogram as the original TH1. For certain type of
00239   //                                 bin contents (containing efficiencies, asymmetries, or ratio is general)
00240   //                                 you should import the absolute value and set impDens to kFALSE
00241   //                                 
00242   //
00243   // Weight(Double_t)          -- Apply given weight factor when importing histograms
00244   //
00245   // Index(RooCategory&)       -- Prepare import of multiple TH1/1/2/3 into a N+1 dimensional RooDataHist
00246   //                              where the extra discrete dimension labels the source of the imported histogram
00247   //                              If the index category defines states for which no histogram is be imported
00248   //                              the corresponding bins will be left empty.
00249   //                              
00250   // Import(const char*, TH1&) -- Import a THx to be associated with the given state name of the index category
00251   //                              specified in Index(). If the given state name is not yet defined in the index
00252   //                              category it will be added on the fly. The import command can be specified
00253   //                              multiple times. 
00254   // Import(map<string,TH1*>&) -- As above, but allows specification of many imports in a single operation
00255   //                              
00256 
00257   // Initialize datastore
00258   _dstore = new RooTreeDataStore(name,title,_vars) ; 
00259 
00260   // Define configuration for this method
00261   RooCmdConfig pc(Form("RooDataHist::ctor(%s)",GetName())) ;
00262   pc.defineObject("impHist","ImportHisto",0) ;
00263   pc.defineInt("impDens","ImportHisto",0) ;
00264   pc.defineObject("indexCat","IndexCat",0) ;
00265   pc.defineObject("impSliceHist","ImportHistoSlice",0,0,kTRUE) ; // array
00266   pc.defineString("impSliceState","ImportHistoSlice",0,"",kTRUE) ; // array
00267   pc.defineObject("impSliceDHist","ImportDataHistSlice",0,0,kTRUE) ; // array
00268   pc.defineString("impSliceDState","ImportDataHistSlice",0,"",kTRUE) ; // array
00269   pc.defineDouble("weight","Weight",0,1) ; 
00270   pc.defineObject("dummy1","ImportDataHistSliceMany",0) ;
00271   pc.defineObject("dummy2","ImportHistoSliceMany",0) ;
00272   pc.defineMutex("ImportHisto","ImportHistoSlice","ImportDataHistSlice") ;
00273   pc.defineDependency("ImportHistoSlice","IndexCat") ;
00274   pc.defineDependency("ImportDataHistSlice","IndexCat") ;
00275 
00276   RooLinkedList l ;
00277   l.Add((TObject*)&arg1) ;  l.Add((TObject*)&arg2) ;  
00278   l.Add((TObject*)&arg3) ;  l.Add((TObject*)&arg4) ;
00279   l.Add((TObject*)&arg5) ;  l.Add((TObject*)&arg6) ;  
00280   l.Add((TObject*)&arg7) ;  l.Add((TObject*)&arg8) ;
00281 
00282   // Process & check varargs 
00283   pc.process(l) ;
00284   if (!pc.ok(kTRUE)) {
00285     assert(0) ;
00286     return ;
00287   }
00288 
00289   TH1* impHist = static_cast<TH1*>(pc.getObject("impHist")) ;
00290   Bool_t impDens = pc.getInt("impDens") ;
00291   Double_t initWgt = pc.getDouble("weight") ;
00292   const char* impSliceNames = pc.getString("impSliceState","",kTRUE) ;
00293   const RooLinkedList& impSliceHistos = pc.getObjectList("impSliceHist") ;
00294   RooCategory* indexCat = static_cast<RooCategory*>(pc.getObject("indexCat")) ;
00295   const char* impSliceDNames = pc.getString("impSliceDState","",kTRUE) ;
00296   const RooLinkedList& impSliceDHistos = pc.getObjectList("impSliceDHist") ;
00297 
00298 
00299   if (impHist) {
00300     
00301     // Initialize importing contents from TH1
00302     importTH1(vars,*impHist,initWgt, impDens) ;
00303 
00304   } else if (indexCat) {
00305 
00306 
00307     if (impSliceHistos.GetSize()>0) {
00308 
00309       // Initialize importing mapped set of TH1s
00310       map<string,TH1*> hmap ;
00311       char tmp[1024] ;
00312       strlcpy(tmp,impSliceNames,1024) ;
00313       char* token = strtok(tmp,",") ;
00314       TIterator* hiter = impSliceHistos.MakeIterator() ;
00315       while(token) {
00316         hmap[token] = (TH1*) hiter->Next() ;
00317         token = strtok(0,",") ;
00318       }
00319       importTH1Set(vars,*indexCat,hmap,initWgt,kTRUE) ;
00320     } else {
00321 
00322       // Initialize importing mapped set of RooDataHists
00323       map<string,RooDataHist*> dmap ;
00324       char tmp[1024] ;
00325       strlcpy(tmp,impSliceDNames,1024) ;
00326       char* token = strtok(tmp,",") ;
00327       TIterator* hiter = impSliceDHistos.MakeIterator() ;
00328       while(token) {
00329         dmap[token] = (RooDataHist*) hiter->Next() ;
00330         token = strtok(0,",") ;
00331       }
00332       importDHistSet(vars,*indexCat,dmap,initWgt) ;
00333 
00334     }
00335 
00336 
00337   } else {
00338 
00339     // Initialize empty
00340     initialize() ;
00341     appendToDir(this,kTRUE) ;    
00342 
00343   }
00344 
00345   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00346 
00347 }
00348 
00349 
00350 
00351 
00352 //_____________________________________________________________________________
00353 void RooDataHist::importTH1(const RooArgList& vars, TH1& histo, Double_t wgt, Bool_t doDensityCorrection) 
00354 {
00355   // Import data from given TH1/2/3 into this RooDataHist
00356 
00357   // Adjust binning of internal observables to match that of input THx
00358   Int_t offset[3] ;
00359   adjustBinning(vars,histo,offset) ;
00360   
00361   // Initialize internal data structure
00362   initialize() ;
00363   appendToDir(this,kTRUE) ;
00364 
00365   // Define x,y,z as 1st, 2nd and 3rd observable
00366   RooRealVar* xvar = (RooRealVar*) _vars.find(vars.at(0)->GetName()) ;
00367   RooRealVar* yvar = (RooRealVar*) (vars.at(1) ? _vars.find(vars.at(1)->GetName()) : 0 ) ;
00368   RooRealVar* zvar = (RooRealVar*) (vars.at(2) ? _vars.find(vars.at(2)->GetName()) : 0 ) ;
00369 
00370   // Transfer contents
00371   Int_t xmin(0),ymin(0),zmin(0) ;
00372   RooArgSet vset(*xvar) ;
00373   Double_t volume = xvar->getMax()-xvar->getMin() ;
00374   xmin = offset[0] ;
00375   if (yvar) {
00376     vset.add(*yvar) ;
00377     ymin = offset[1] ;
00378     volume *= (yvar->getMax()-yvar->getMin()) ;
00379   }
00380   if (zvar) {
00381     vset.add(*zvar) ;
00382     zmin = offset[2] ;
00383     volume *= (zvar->getMax()-zvar->getMin()) ;
00384   }
00385   Double_t avgBV = volume / numEntries() ;
00386 
00387   Int_t ix(0),iy(0),iz(0) ;
00388   for (ix=0 ; ix < xvar->getBins() ; ix++) {
00389     xvar->setBin(ix) ;
00390     if (yvar) {
00391       for (iy=0 ; iy < yvar->getBins() ; iy++) {
00392         yvar->setBin(iy) ;
00393         if (zvar) {
00394           for (iz=0 ; iz < zvar->getBins() ; iz++) {
00395             zvar->setBin(iz) ;
00396             Double_t bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
00397             add(vset,bv*histo.GetBinContent(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,bv*TMath::Power(histo.GetBinError(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,2)) ;
00398           }
00399         } else {
00400           Double_t bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
00401           add(vset,bv*histo.GetBinContent(ix+1+xmin,iy+1+ymin)*wgt,bv*TMath::Power(histo.GetBinError(ix+1+xmin,iy+1+ymin)*wgt,2)) ;
00402         }
00403       }
00404     } else {
00405       Double_t bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
00406       add(vset,bv*histo.GetBinContent(ix+1+xmin)*wgt,bv*TMath::Power(histo.GetBinError(ix+1+xmin)*wgt,2)) ;         
00407     }
00408   }  
00409   
00410 }
00411 
00412 
00413 
00414 
00415 
00416 //_____________________________________________________________________________
00417 void RooDataHist::importTH1Set(const RooArgList& vars, RooCategory& indexCat, map<string,TH1*> hmap, Double_t wgt, Bool_t doDensityCorrection) 
00418 {
00419   // Import data from given set of TH1/2/3 into this RooDataHist. The category indexCat labels the sources
00420   // in the constructed RooDataHist. The stl map provides the mapping between the indexCat state labels
00421   // and the import source
00422 
00423   RooCategory* icat = (RooCategory*) _vars.find(indexCat.GetName()) ;
00424 
00425   TH1* histo(0) ;  
00426   Bool_t init(kFALSE) ;
00427   for (map<string,TH1*>::iterator hiter = hmap.begin() ; hiter!=hmap.end() ; ++hiter) {
00428     // Store pointer to first histogram from which binning specification will be taken
00429     if (!histo) {
00430       histo = hiter->second ;
00431     }
00432     // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
00433     if (!indexCat.lookupType(hiter->first.c_str())) {
00434       indexCat.defineType(hiter->first.c_str()) ;
00435       coutI(InputArguments) << "RooDataHist::importTH1Set(" << GetName() << ") defining state \"" << hiter->first << "\" in index category " << indexCat.GetName() << endl ;
00436     }
00437     if (!icat->lookupType(hiter->first.c_str())) {      
00438       icat->defineType(hiter->first.c_str()) ;
00439     }
00440   }
00441 
00442   // Check consistency in number of dimensions
00443   if (histo && (vars.getSize() != histo->GetDimension())) {
00444     coutE(InputArguments) << "RooDataHist::ctor(" << GetName() << ") ERROR: dimension of input histogram must match "
00445                           << "number of continuous variables" << endl ;
00446     assert(0) ; 
00447   }
00448   
00449   // Copy bins and ranges from THx to dimension observables
00450   Int_t offset[3] ;
00451   adjustBinning(vars,*histo,offset) ;
00452   
00453   // Initialize internal data structure
00454   if (!init) {
00455     initialize() ;
00456     appendToDir(this,kTRUE) ;
00457     init = kTRUE ;
00458   }
00459 
00460   // Define x,y,z as 1st, 2nd and 3rd observable
00461   RooRealVar* xvar = (RooRealVar*) _vars.find(vars.at(0)->GetName()) ;
00462   RooRealVar* yvar = (RooRealVar*) (vars.at(1) ? _vars.find(vars.at(1)->GetName()) : 0 ) ;
00463   RooRealVar* zvar = (RooRealVar*) (vars.at(2) ? _vars.find(vars.at(2)->GetName()) : 0 ) ;
00464 
00465   // Transfer contents
00466   Int_t xmin(0),ymin(0),zmin(0) ;
00467   RooArgSet vset(*xvar) ;
00468   Double_t volume = xvar->getMax()-xvar->getMin() ;
00469   xmin = offset[0] ;
00470   if (yvar) {
00471     vset.add(*yvar) ;
00472     ymin = offset[1] ;
00473     volume *= (yvar->getMax()-yvar->getMin()) ;
00474   }
00475   if (zvar) {
00476     vset.add(*zvar) ;
00477     zmin = offset[2] ;
00478     volume *= (zvar->getMax()-zvar->getMin()) ;
00479   }
00480   Double_t avgBV = volume / numEntries() ;
00481   
00482   Int_t ic(0),ix(0),iy(0),iz(0) ;
00483   for (ic=0 ; ic < icat->numBins(0) ; ic++) {
00484     icat->setBin(ic) ;
00485     histo = hmap[icat->getLabel()] ;
00486     for (ix=0 ; ix < xvar->getBins() ; ix++) {
00487       xvar->setBin(ix) ;
00488       if (yvar) {
00489         for (iy=0 ; iy < yvar->getBins() ; iy++) {
00490           yvar->setBin(iy) ;
00491           if (zvar) {
00492             for (iz=0 ; iz < zvar->getBins() ; iz++) {
00493               zvar->setBin(iz) ;
00494               Double_t bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
00495               add(vset,bv*histo->GetBinContent(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,bv*TMath::Power(histo->GetBinError(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,2)) ;
00496             }
00497           } else {
00498             Double_t bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
00499             add(vset,bv*histo->GetBinContent(ix+1+xmin,iy+1+ymin)*wgt,bv*TMath::Power(histo->GetBinError(ix+1+xmin,iy+1+ymin)*wgt,2)) ;
00500           }
00501         }
00502       } else {
00503         Double_t bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
00504         add(vset,bv*histo->GetBinContent(ix+1+xmin)*wgt,bv*TMath::Power(histo->GetBinError(ix+1+xmin)*wgt,2)) ;     
00505       }
00506     }  
00507   }
00508 
00509 }
00510 
00511 
00512 
00513 //_____________________________________________________________________________
00514 void RooDataHist::importDHistSet(const RooArgList& /*vars*/, RooCategory& indexCat, std::map<std::string,RooDataHist*> dmap, Double_t initWgt) 
00515 {
00516   // Import data from given set of TH1/2/3 into this RooDataHist. The category indexCat labels the sources
00517   // in the constructed RooDataHist. The stl map provides the mapping between the indexCat state labels
00518   // and the import source
00519 
00520   RooCategory* icat = (RooCategory*) _vars.find(indexCat.GetName()) ;
00521 
00522   for (map<string,RooDataHist*>::iterator diter = dmap.begin() ; diter!=dmap.end() ; ++diter) {
00523 
00524     // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
00525     if (!indexCat.lookupType(diter->first.c_str())) {
00526       indexCat.defineType(diter->first.c_str()) ;
00527       coutI(InputArguments) << "RooDataHist::importDHistSet(" << GetName() << ") defining state \"" << diter->first << "\" in index category " << indexCat.GetName() << endl ;
00528     }
00529     if (!icat->lookupType(diter->first.c_str())) {      
00530       icat->defineType(diter->first.c_str()) ;
00531     }
00532   }
00533 
00534   initialize() ;
00535   appendToDir(this,kTRUE) ;  
00536 
00537 
00538   for (map<string,RooDataHist*>::iterator diter = dmap.begin() ; diter!=dmap.end() ; ++diter) {
00539 
00540     RooDataHist* dhist = diter->second ;
00541 
00542     icat->setLabel(diter->first.c_str()) ;
00543 
00544     // Transfer contents
00545     for (Int_t i=0 ; i<dhist->numEntries() ; i++) {
00546       _vars = *dhist->get(i) ;
00547       add(_vars,dhist->weight()*initWgt, pow(dhist->weightError(SumW2),2) ) ;
00548     }
00549 
00550   }
00551 }
00552 
00553 //_____________________________________________________________________________
00554 void RooDataHist::adjustBinning(const RooArgList& vars, TH1& href, Int_t* offset) 
00555 {
00556   // Adjust binning specification on first and optionally second and third
00557   // observable to binning in given reference TH1. Used by constructors
00558   // that import data from an external TH1
00559 
00560   // X
00561   RooRealVar* xvar = (RooRealVar*) _vars.find(vars.at(0)->GetName()) ;
00562   if (!dynamic_cast<RooRealVar*>(xvar)) {
00563     coutE(InputArguments) << "RooDataHist::adjustBinning(" << GetName() << ") ERROR: dimension " << xvar->GetName() << " must be real" << endl ;
00564     assert(0) ;
00565   }
00566 
00567   Double_t xlo = ((RooRealVar*)vars.at(0))->getMin() ;
00568   Double_t xhi = ((RooRealVar*)vars.at(0))->getMax() ;
00569   Int_t xmin(0) ;
00570   if (href.GetXaxis()->GetXbins()->GetArray()) {
00571 
00572     RooBinning xbins(href.GetNbinsX(),href.GetXaxis()->GetXbins()->GetArray()) ;
00573 
00574     // Adjust xlo/xhi to nearest boundary
00575     Double_t xloAdj = xbins.binLow(xbins.binNumber(xlo+1e-6)) ;
00576     Double_t xhiAdj = xbins.binHigh(xbins.binNumber(xhi-1e-6)) ;
00577     xbins.setRange(xloAdj,xhiAdj) ;
00578     ((RooRealVar*)vars.at(0))->setBinning(xbins) ;
00579     if (fabs(xloAdj-xlo)>1e-6||fabs(xhiAdj-xhi)) {
00580       coutI(DataHandling) << "RooDataHist::adjustBinning(" << GetName() << "): fit range of variable " << xvar->GetName() << " expanded to nearest bin boundaries: [" 
00581                           << xlo << "," << xhi << "] --> [" << xloAdj << "," << xhiAdj << "]" << endl ;
00582     }
00583 
00584     xvar->setBinning(xbins) ;
00585     xmin = xbins.rawBinNumber(xloAdj+1e-6) ;
00586     if (offset) {
00587       offset[0] = xmin ;
00588     }
00589 
00590   } else {
00591     RooBinning xbins(href.GetXaxis()->GetXmin(),href.GetXaxis()->GetXmax()) ;
00592     xbins.addUniform(href.GetNbinsX(),href.GetXaxis()->GetXmin(),href.GetXaxis()->GetXmax()) ;
00593 
00594     // Adjust xlo/xhi to nearest boundary
00595     Double_t xloAdj = xbins.binLow(xbins.binNumber(xlo+1e-6)) ;
00596     Double_t xhiAdj = xbins.binHigh(xbins.binNumber(xhi-1e-6)) ;
00597     xbins.setRange(xloAdj,xhiAdj) ;
00598     ((RooRealVar*)vars.at(0))->setRange(xloAdj,xhiAdj) ;
00599     if (fabs(xloAdj-xlo)>1e-6||fabs(xhiAdj-xhi)) {
00600       coutI(DataHandling) << "RooDataHist::adjustBinning(" << GetName() << "): fit range of variable " << xvar->GetName() << " expanded to nearest bin boundaries: [" 
00601                           << xlo << "," << xhi << "] --> [" << xloAdj << "," << xhiAdj << "]" << endl ;
00602     }
00603 
00604 
00605     RooUniformBinning xbins2(xloAdj,xhiAdj,xbins.numBins()) ;
00606     xvar->setBinning(xbins2) ;
00607     xmin = xbins.rawBinNumber(xloAdj+1e-6) ;
00608     if (offset) {
00609       offset[0] = xmin ;
00610     }
00611   }
00612 
00613 
00614 
00615   // Y
00616   RooRealVar* yvar = (RooRealVar*) (vars.at(1) ? _vars.find(vars.at(1)->GetName()) : 0 ) ;
00617   Int_t ymin(0) ;
00618   if (yvar) {
00619     Double_t ylo = ((RooRealVar*)vars.at(1))->getMin() ;
00620     Double_t yhi = ((RooRealVar*)vars.at(1))->getMax() ;
00621 
00622     if (!dynamic_cast<RooRealVar*>(yvar)) {
00623       coutE(InputArguments) << "RooDataHist::adjustBinning(" << GetName() << ") ERROR: dimension " << yvar->GetName() << " must be real" << endl ;
00624       assert(0) ;
00625     }
00626 
00627     if (href.GetYaxis()->GetXbins()->GetArray()) {
00628 
00629       RooBinning ybins(href.GetNbinsY(),href.GetYaxis()->GetXbins()->GetArray()) ;
00630       
00631       // Adjust ylo/yhi to nearest boundary
00632       Double_t yloAdj = ybins.binLow(ybins.binNumber(ylo+1e-6)) ;
00633       Double_t yhiAdj = ybins.binHigh(ybins.binNumber(yhi-1e-6)) ;
00634       ybins.setRange(yloAdj,yhiAdj) ;
00635       ((RooRealVar*)vars.at(1))->setBinning(ybins) ;
00636       if (fabs(yloAdj-ylo)>1e-6||fabs(yhiAdj-yhi)) {
00637         coutI(DataHandling) << "RooDataHist::adjustBinning(" << GetName() << "): fit range of variable " << yvar->GetName() << " expanded to nearest bin boundaries: [" 
00638                             << ylo << "," << yhi << "] --> [" << yloAdj << "," << yhiAdj << "]" << endl ;
00639       }
00640 
00641       yvar->setBinning(ybins) ;
00642       ymin = ybins.rawBinNumber(yloAdj+1e-6) ;
00643       if (offset) {
00644         offset[1] = ymin ;
00645       }
00646 
00647     } else {
00648 
00649       RooBinning ybins(href.GetYaxis()->GetXmin(),href.GetYaxis()->GetXmax()) ;
00650       ybins.addUniform(href.GetNbinsY(),href.GetYaxis()->GetXmin(),href.GetYaxis()->GetXmax()) ;
00651       
00652       // Adjust ylo/yhi to nearest boundary
00653       Double_t yloAdj = ybins.binLow(ybins.binNumber(ylo+1e-6)) ;
00654       Double_t yhiAdj = ybins.binHigh(ybins.binNumber(yhi-1e-6)) ;
00655       ybins.setRange(yloAdj,yhiAdj) ;
00656       ((RooRealVar*)vars.at(1))->setRange(yloAdj,yhiAdj) ;
00657       if (fabs(yloAdj-ylo)>1e-6||fabs(yhiAdj-yhi)) {
00658         coutI(DataHandling) << "RooDataHist::adjustBinning(" << GetName() << "): fit range of variable " << yvar->GetName() << " expanded to nearest bin boundaries: [" 
00659                             << ylo << "," << yhi << "] --> [" << yloAdj << "," << yhiAdj << "]" << endl ;
00660       }
00661       
00662       RooUniformBinning ybins2(yloAdj,yhiAdj,ybins.numBins()) ;
00663       yvar->setBinning(ybins2) ;
00664       ymin = ybins.rawBinNumber(yloAdj+1e-6) ;
00665       if (offset) {
00666         offset[1] = ymin ;
00667       }
00668 
00669     }    
00670   }
00671   
00672   // Z
00673   RooRealVar* zvar = (RooRealVar*) (vars.at(2) ? _vars.find(vars.at(2)->GetName()) : 0 ) ;
00674   Int_t zmin(0) ;
00675   if (zvar) {
00676     Double_t zlo = ((RooRealVar*)vars.at(2))->getMin() ;
00677     Double_t zhi = ((RooRealVar*)vars.at(2))->getMax() ;
00678 
00679     if (!dynamic_cast<RooRealVar*>(zvar)) {
00680       coutE(InputArguments) << "RooDataHist::adjustBinning(" << GetName() << ") ERROR: dimension " << zvar->GetName() << " must be real" << endl ;
00681       assert(0) ;
00682     }
00683 
00684     if (href.GetZaxis()->GetXbins()->GetArray()) {
00685 
00686       RooBinning zbins(href.GetNbinsZ(),href.GetZaxis()->GetXbins()->GetArray()) ;
00687       
00688       // Adjust zlo/zhi to nearest boundary
00689       Double_t zloAdj = zbins.binLow(zbins.binNumber(zlo+1e-6)) ;
00690       Double_t zhiAdj = zbins.binHigh(zbins.binNumber(zhi-1e-6)) ;
00691       zbins.setRange(zloAdj,zhiAdj) ;
00692       ((RooRealVar*)vars.at(2))->setBinning(zbins) ;
00693       if (fabs(zloAdj-zlo)>1e-6||fabs(zhiAdj-zhi)) {
00694         coutI(DataHandling) << "RooDataHist::adjustBinning(" << GetName() << "): fit range of variable " << zvar->GetName() << " expanded to nearest bin boundaries: [" 
00695                             << zlo << "," << zhi << "] --> [" << zloAdj << "," << zhiAdj << "]" << endl ;
00696       }
00697       
00698       zvar->setBinning(zbins) ;
00699       zmin = zbins.rawBinNumber(zloAdj+1e-6) ;
00700       if (offset) {
00701         offset[2] = zmin ;
00702       }
00703       
00704     } else {
00705 
00706       RooBinning zbins(href.GetZaxis()->GetXmin(),href.GetZaxis()->GetXmax()) ;
00707       zbins.addUniform(href.GetNbinsZ(),href.GetZaxis()->GetXmin(),href.GetZaxis()->GetXmax()) ;
00708       
00709       // Adjust zlo/zhi to nearest boundary
00710       Double_t zloAdj = zbins.binLow(zbins.binNumber(zlo+1e-6)) ;
00711       Double_t zhiAdj = zbins.binHigh(zbins.binNumber(zhi-1e-6)) ;
00712       zbins.setRange(zloAdj,zhiAdj) ;
00713       ((RooRealVar*)vars.at(2))->setRange(zloAdj,zhiAdj) ;
00714       if (fabs(zloAdj-zlo)>1e-6||fabs(zhiAdj-zhi)) {
00715         coutI(DataHandling) << "RooDataHist::adjustBinning(" << GetName() << "): fit range of variable " << zvar->GetName() << " expanded to nearest bin boundaries: [" 
00716                             << zlo << "," << zhi << "] --> [" << zloAdj << "," << zhiAdj << "]" << endl ;
00717       }
00718       
00719       RooUniformBinning zbins2(zloAdj,zhiAdj,zbins.numBins()) ;
00720       zvar->setBinning(zbins2) ;
00721       zmin = zbins.rawBinNumber(zloAdj+1e-6) ;
00722       if (offset) {
00723         offset[2] = zmin ;
00724       }
00725     }
00726   }
00727 
00728 }
00729 
00730 
00731 
00732 
00733 
00734 //_____________________________________________________________________________
00735 void RooDataHist::initialize(const char* binningName, Bool_t fillTree)
00736 {
00737   // Initialization procedure: allocate weights array, calculate
00738   // multipliers needed for N-space to 1-dim array jump table,
00739   // and fill the internal tree with all bin center coordinates
00740 
00741 
00742   // Save real dimensions of dataset separately
00743   RooAbsArg* real ;
00744   _iterator->Reset() ;
00745   while((real=(RooAbsArg*)_iterator->Next())) {
00746     if (dynamic_cast<RooAbsReal*>(real)) _realVars.add(*real) ;
00747   }
00748   _realIter = _realVars.createIterator() ;
00749 
00750   // Fill array of LValue pointers to variables
00751   _iterator->Reset() ;
00752   RooAbsArg* rvarg ;
00753   while((rvarg=(RooAbsArg*)_iterator->Next())) {
00754     if (binningName) {
00755       RooRealVar* rrv = dynamic_cast<RooRealVar*>(rvarg) ; 
00756       if (rrv) {
00757         rrv->setBinning(rrv->getBinning(binningName)) ;
00758       }
00759     }
00760     // coverity[FORWARD_NULL]
00761     _lvvars.push_back(dynamic_cast<RooAbsLValue*>(rvarg)) ;    
00762     // coverity[FORWARD_NULL]
00763     const RooAbsBinning* binning = dynamic_cast<RooAbsLValue*>(rvarg)->getBinningPtr(0) ;
00764     _lvbins.push_back(binning ? binning->clone() : 0) ;    
00765   }
00766 
00767   
00768   // Allocate coefficients array
00769   _idxMult.resize(_vars.getSize()) ;
00770 
00771   _arrSize = 1 ;
00772   _iterator->Reset() ;
00773   RooAbsLValue* arg ;
00774   Int_t n(0), i ;
00775   while((arg=dynamic_cast<RooAbsLValue*>(_iterator->Next()))) {
00776     
00777     // Calculate sub-index multipliers for master index
00778     for (i=0 ; i<n ; i++) {
00779       _idxMult[i] *= arg->numBins() ;
00780     }
00781     _idxMult[n++] = 1 ;
00782 
00783     // Calculate dimension of weight array
00784     _arrSize *= arg->numBins() ;
00785   }  
00786 
00787   // Allocate and initialize weight array if necessary
00788   if (!_wgt) {
00789     _wgt = new Double_t[_arrSize] ;
00790     _errLo = new Double_t[_arrSize] ;
00791     _errHi = new Double_t[_arrSize] ;
00792     _sumw2 = new Double_t[_arrSize] ;
00793     _binv = new Double_t[_arrSize] ;
00794 
00795     // Refill array pointers in data store when reading
00796     // from Streamer
00797     if (fillTree==kFALSE) {
00798       ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;      
00799     }
00800     
00801     for (i=0 ; i<_arrSize ; i++) {
00802       _wgt[i] = 0 ;
00803       _errLo[i] = -1 ;
00804       _errHi[i] = -1 ;
00805       _sumw2[i] = 0 ;
00806     }
00807   }
00808 
00809   if (!fillTree) return ;
00810 
00811   // Fill TTree with bin center coordinates
00812   // Calculate plot bins of components from master index
00813 
00814   Int_t ibin ;
00815   for (ibin=0 ; ibin<_arrSize ; ibin++) {
00816     _iterator->Reset() ;
00817     RooAbsLValue* arg2 ;
00818     Int_t j(0), idx(0), tmp(ibin) ;
00819     Double_t theBinVolume(1) ;
00820     while((arg2=dynamic_cast<RooAbsLValue*>(_iterator->Next()))) {
00821       idx  = tmp / _idxMult[j] ;
00822       tmp -= idx*_idxMult[j++] ;
00823       RooAbsLValue* arglv = dynamic_cast<RooAbsLValue*>(arg2) ;
00824       arglv->setBin(idx) ;
00825       theBinVolume *= arglv->getBinWidth(idx) ;
00826     }
00827     _binv[ibin] = theBinVolume ;
00828     fill() ;
00829   }
00830 
00831 
00832 }
00833 
00834 
00835 
00836 //_____________________________________________________________________________
00837 RooDataHist::RooDataHist(const RooDataHist& other, const char* newname) :
00838   RooAbsData(other,newname), RooDirItem(), _idxMult(other._idxMult), _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(other._pbinvCacheMgr,0)
00839 {
00840   // Copy constructor
00841 
00842   Int_t i ;
00843 
00844   // Allocate and initialize weight array 
00845   _arrSize = other._arrSize ;
00846   _wgt = new Double_t[_arrSize] ;
00847   _errLo = new Double_t[_arrSize] ;
00848   _errHi = new Double_t[_arrSize] ;
00849   _binv = new Double_t[_arrSize] ;
00850   _sumw2 = new Double_t[_arrSize] ;
00851   for (i=0 ; i<_arrSize ; i++) {
00852     _wgt[i] = other._wgt[i] ;
00853     _errLo[i] = other._errLo[i] ;
00854     _errHi[i] = other._errHi[i] ;
00855     _sumw2[i] = other._sumw2[i] ;
00856     _binv[i] = other._binv[i] ;
00857   }  
00858 
00859   // Save real dimensions of dataset separately
00860   RooAbsArg* arg ;
00861   _iterator->Reset() ;
00862   while((arg=(RooAbsArg*)_iterator->Next())) {
00863     if (dynamic_cast<RooAbsReal*>(arg)) _realVars.add(*arg) ;
00864   }
00865   _realIter = _realVars.createIterator() ;
00866 
00867   // Fill array of LValue pointers to variables
00868   _iterator->Reset() ;
00869   RooAbsArg* rvarg ;
00870   while((rvarg=(RooAbsArg*)_iterator->Next())) {
00871     // coverity[FORWARD_NULL]
00872     _lvvars.push_back(dynamic_cast<RooAbsLValue*>(rvarg)) ;
00873     // coverity[FORWARD_NULL]
00874     const RooAbsBinning* binning = dynamic_cast<RooAbsLValue*>(rvarg)->getBinningPtr(0) ;
00875     _lvbins.push_back(binning ? binning->clone() : 0) ;    
00876   }
00877 
00878   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00879 
00880  appendToDir(this,kTRUE) ;
00881 }
00882 
00883 
00884 
00885 //_____________________________________________________________________________
00886 RooDataHist::RooDataHist(const char* name, const char* title, RooDataHist* h, const RooArgSet& varSubset, 
00887                          const RooFormulaVar* cutVar, const char* cutRange, Int_t nStart, Int_t nStop, Bool_t copyCache) :
00888   RooAbsData(name,title,varSubset),
00889   _binValid(0), _curWeight(0), _curVolume(1), _pbinv(0), _pbinvCacheMgr(0,10)
00890 {
00891   // Constructor of a data hist from (part of) an existing data hist. The dimensions
00892   // of the data set are defined by the 'vars' RooArgSet, which can be identical
00893   // to 'dset' dimensions, or a subset thereof. Reduced dimensions will be projected
00894   // in the output data hist. The optional 'cutVar' formula variable can used to 
00895   // select the subset of bins to be copied.
00896   //
00897   // For most uses the RooAbsData::reduce() wrapper function, which uses this constructor, 
00898   // is the most convenient way to create a subset of an existing data  
00899 
00900   // Initialize datastore
00901   _dstore = new RooTreeDataStore(name,title,*h->_dstore,_vars,cutVar,cutRange,nStart,nStop,copyCache) ;
00902   
00903   initialize(0,kFALSE) ;
00904 
00905   ((RooTreeDataStore*)_dstore)->setExternalWeightArray(_wgt,_errLo,_errHi,_sumw2) ;
00906 
00907   // Copy weight array etc
00908   Int_t i ;
00909   for (i=0 ; i<_arrSize ; i++) {
00910     _wgt[i] = h->_wgt[i] ;
00911     _errLo[i] = h->_errLo[i] ;
00912     _errHi[i] = h->_errHi[i] ;
00913     _sumw2[i] = h->_sumw2[i] ;
00914     _binv[i] = h->_binv[i] ;
00915   }  
00916 
00917 
00918   appendToDir(this,kTRUE) ;
00919 }
00920 
00921 
00922 //_____________________________________________________________________________
00923 RooAbsData* RooDataHist::cacheClone(const RooAbsArg* newCacheOwner, const RooArgSet* newCacheVars, const char* newName) 
00924 {
00925   // Construct a clone of this dataset that contains only the cached variables
00926   checkInit() ;
00927 
00928   RooDataHist* dhist = new RooDataHist(newName?newName:GetName(),GetTitle(),this,*get(),0,0,0,2000000000,kTRUE) ; 
00929 
00930   RooArgSet* selCacheVars = (RooArgSet*) newCacheVars->selectCommon(dhist->_cachedVars) ;
00931   dhist->attachCache(newCacheOwner, *selCacheVars) ;
00932   delete selCacheVars ;
00933 
00934   return dhist ;
00935 }
00936 
00937 
00938 
00939 //_____________________________________________________________________________
00940 RooAbsData* RooDataHist::reduceEng(const RooArgSet& varSubset, const RooFormulaVar* cutVar, const char* cutRange, 
00941                                    Int_t nStart, Int_t nStop, Bool_t /*copyCache*/)
00942 {
00943   // Implementation of RooAbsData virtual method that drives the RooAbsData::reduce() methods
00944 
00945   checkInit() ;
00946   RooArgSet* myVarSubset = (RooArgSet*) _vars.selectCommon(varSubset) ;
00947   RooDataHist *rdh = new RooDataHist(GetName(), GetTitle(), *myVarSubset) ;
00948   delete myVarSubset ;
00949 
00950   RooFormulaVar* cloneVar = 0;
00951   RooArgSet* tmp(0) ;
00952   if (cutVar) {
00953     // Deep clone cutVar and attach clone to this dataset
00954     tmp = (RooArgSet*) RooArgSet(*cutVar).snapshot() ;
00955     if (!tmp) {
00956       coutE(DataHandling) << "RooDataHist::reduceEng(" << GetName() << ") Couldn't deep-clone cut variable, abort," << endl ;
00957       return 0 ;
00958     }
00959     cloneVar = (RooFormulaVar*) tmp->find(cutVar->GetName()) ;
00960     cloneVar->attachDataSet(*this) ;
00961   }
00962 
00963   Int_t i ;
00964   Double_t lo,hi ;
00965   Int_t nevt = nStop < numEntries() ? nStop : numEntries() ;
00966   TIterator* vIter = get()->createIterator() ;
00967   for (i=nStart ; i<nevt ; i++) {
00968     const RooArgSet* row = get(i) ;
00969 
00970     Bool_t doSelect(kTRUE) ;
00971     if (cutRange) {
00972       RooAbsArg* arg ;
00973       vIter->Reset() ;
00974       while((arg=(RooAbsArg*)vIter->Next())) {  
00975         if (!arg->inRange(cutRange)) {
00976           doSelect = kFALSE ;
00977           break ;
00978         }
00979       }
00980     }
00981     if (!doSelect) continue ;
00982 
00983     if (!cloneVar || cloneVar->getVal()) {
00984       weightError(lo,hi,SumW2) ;
00985       rdh->add(*row,weight(),lo*lo) ;
00986     }
00987   }
00988   delete vIter ;
00989 
00990   if (cloneVar) {
00991     delete tmp ;
00992   } 
00993   
00994     return rdh ;
00995   }
00996 
00997 
00998 
00999 //_____________________________________________________________________________
01000 RooDataHist::~RooDataHist() 
01001 {
01002   // Destructor
01003 
01004   if (_wgt) delete[] _wgt ;
01005   if (_errLo) delete[] _errLo ;
01006   if (_errHi) delete[] _errHi ;
01007   if (_sumw2) delete[] _sumw2 ;
01008   if (_binv) delete[] _binv ;
01009   if (_realIter) delete _realIter ;
01010   if (_binValid) delete[] _binValid ;
01011   list<const RooAbsBinning*>::iterator iter = _lvbins.begin() ;
01012   while(iter!=_lvbins.end()) {
01013     delete *iter ;
01014     iter++ ;
01015   }
01016 
01017    removeFromDir(this) ;
01018 }
01019 
01020 
01021 
01022 
01023 //_____________________________________________________________________________
01024 Int_t RooDataHist::getIndex(const RooArgSet& coord)
01025 {
01026   checkInit() ;
01027   _vars = coord ;
01028   return calcTreeIndex() ;
01029 }
01030 
01031 
01032 
01033 
01034 //_____________________________________________________________________________
01035 Int_t RooDataHist::calcTreeIndex() const 
01036 {
01037   // Calculate the index for the weights array corresponding to 
01038   // to the bin enclosing the current coordinates of the internal argset
01039 
01040   Int_t masterIdx(0), i(0) ;
01041   list<RooAbsLValue*>::const_iterator iter = _lvvars.begin() ;
01042   list<const RooAbsBinning*>::const_iterator biter = _lvbins.begin() ;
01043   for (;iter!=_lvvars.end() ; ++iter) {
01044     const RooAbsBinning* binning = (*biter) ;
01045     masterIdx += _idxMult[i++]*(*iter)->getBin(binning) ;
01046     biter++ ;
01047   }
01048   return masterIdx ;
01049 }
01050 
01051 
01052 
01053 //_____________________________________________________________________________
01054 void RooDataHist::dump2() 
01055 {  
01056   // Debug stuff, should go...
01057   cout << "_arrSize = " << _arrSize << endl ;
01058   for (Int_t i=0 ; i<_arrSize ; i++) {
01059     cout << "wgt[" << i << "] = " << _wgt[i] << "sumw2[" << i << "] = " << _sumw2[i] << " vol[" << i << "] = " << _binv[i] << endl ;
01060   }
01061 }
01062 
01063 
01064 
01065 //_____________________________________________________________________________
01066 RooPlot *RooDataHist::plotOn(RooPlot *frame, PlotOpt o) const 
01067 {
01068   // Back end function to plotting functionality. Plot RooDataHist on given
01069   // frame in mode specified by plot options 'o'. The main purpose of
01070   // this function is to match the specified binning on 'o' to the
01071   // internal binning of the plot observable in this RooDataHist.
01072   checkInit() ;
01073   if (o.bins) return RooAbsData::plotOn(frame,o) ;
01074 
01075   if(0 == frame) {
01076     coutE(InputArguments) << ClassName() << "::" << GetName() << ":plotOn: frame is null" << endl;
01077     return 0;
01078   }
01079   RooAbsRealLValue *var= (RooAbsRealLValue*) frame->getPlotVar();
01080   if(0 == var) {
01081     coutE(InputArguments) << ClassName() << "::" << GetName()
01082          << ":plotOn: frame does not specify a plot variable" << endl;
01083     return 0;
01084   }
01085 
01086   RooRealVar* dataVar = (RooRealVar*) _vars.find(var->GetName()) ;
01087   if (!dataVar) {
01088     coutE(InputArguments) << ClassName() << "::" << GetName()
01089          << ":plotOn: dataset doesn't contain plot frame variable" << endl;
01090     return 0;
01091   }
01092 
01093   o.bins = &dataVar->getBinning() ;
01094   o.correctForBinWidth = kFALSE ;
01095   return RooAbsData::plotOn(frame,o) ;
01096 }
01097 
01098 
01099 
01100 
01101 //_____________________________________________________________________________
01102 Double_t RooDataHist::weight(const RooArgSet& bin, Int_t intOrder, Bool_t correctForBinSize, Bool_t cdfBoundaries) 
01103 {
01104   // Return the weight at given coordinates with optional
01105   // interpolation. If intOrder is zero, the weight
01106   // for the bin enclosing the coordinates
01107   // contained in 'bin' is returned. For higher values,
01108   // the result is interpolated in the real dimensions 
01109   // of the dataset
01110 
01111   checkInit() ;
01112 
01113   // Handle illegal intOrder values
01114   if (intOrder<0) {
01115     coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") ERROR: interpolation order must be positive" << endl ;
01116     return 0 ;
01117   }
01118 
01119   // Handle no-interpolation case
01120   if (intOrder==0) {
01121     _vars.assignValueOnly(bin) ;
01122     Int_t idx = calcTreeIndex() ;
01123     //cout << "intOrder 0, idx = " << idx << endl ;
01124     if (correctForBinSize) {
01125       //calculatePartialBinVolume(*get()) ;
01126       return _wgt[idx] / _binv[idx] ;
01127     } else {
01128       return _wgt[idx] ;
01129     }
01130   }
01131 
01132   // Handle all interpolation cases
01133   _vars.assignValueOnly(bin) ;
01134 
01135   Double_t wInt(0) ;
01136   if (_realVars.getSize()==1) {
01137 
01138     // 1-dimensional interpolation
01139     _realIter->Reset() ;
01140     RooRealVar* real=(RooRealVar*)_realIter->Next() ;
01141     const RooAbsBinning* binning = real->getBinningPtr(0) ;
01142     wInt = interpolateDim(*real,binning,((RooAbsReal*)bin.find(real->GetName()))->getVal(), intOrder, correctForBinSize, cdfBoundaries) ;
01143 
01144   } else if (_realVars.getSize()==2) {
01145 
01146     // 2-dimensional interpolation
01147     _realIter->Reset() ;
01148     RooRealVar* realX=(RooRealVar*)_realIter->Next() ;
01149     RooRealVar* realY=(RooRealVar*)_realIter->Next() ;
01150     Double_t xval = ((RooAbsReal*)bin.find(realX->GetName()))->getVal() ;
01151     Double_t yval = ((RooAbsReal*)bin.find(realY->GetName()))->getVal() ;
01152     
01153     Int_t ybinC = realY->getBin() ;
01154     Int_t ybinLo = ybinC-intOrder/2 - ((yval<realY->getBinning().binCenter(ybinC))?1:0) ;
01155     Int_t ybinM = realY->numBins() ;
01156     
01157     Int_t i ;
01158     Double_t yarr[10] ;
01159     Double_t xarr[10] ;
01160     const RooAbsBinning* binning = realX->getBinningPtr(0) ;
01161     for (i=ybinLo ; i<=intOrder+ybinLo ; i++) {
01162       Int_t ibin ;
01163       if (i>=0 && i<ybinM) {
01164         // In range
01165         ibin = i ;
01166         realY->setBin(ibin) ;
01167         xarr[i-ybinLo] = realY->getVal() ;
01168       } else if (i>=ybinM) {
01169         // Overflow: mirror
01170         ibin = 2*ybinM-i-1 ;
01171         realY->setBin(ibin) ;
01172         xarr[i-ybinLo] = 2*realY->getMax()-realY->getVal() ;
01173       } else {
01174         // Underflow: mirror
01175         ibin = -i ;
01176         realY->setBin(ibin) ;
01177         xarr[i-ybinLo] = 2*realY->getMin()-realY->getVal() ;
01178       }
01179       yarr[i-ybinLo] = interpolateDim(*realX,binning,xval,intOrder,correctForBinSize,kFALSE) ;  
01180     }
01181 
01182     if (gDebug>7) {
01183       cout << "RooDataHist interpolating data is" << endl ;
01184       cout << "xarr = " ;
01185       for (int q=0; q<=intOrder ; q++) cout << xarr[q] << " " ;
01186       cout << " yarr = " ;
01187       for (int q=0; q<=intOrder ; q++) cout << yarr[q] << " " ;
01188       cout << endl ;
01189     }
01190     wInt = RooMath::interpolate(xarr,yarr,intOrder+1,yval) ;
01191     
01192   } else {
01193 
01194     // Higher dimensional scenarios not yet implemented
01195     coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") interpolation in " 
01196          << _realVars.getSize() << " dimensions not yet implemented" << endl ;
01197     return weight(bin,0) ;
01198 
01199   }
01200 
01201   // Cut off negative values
01202 //   if (wInt<=0) {
01203 //     wInt=0 ; 
01204 //   }
01205 
01206   //cout << "RooDataHist wInt = " << wInt << endl ;
01207   return wInt ;
01208 }
01209 
01210 
01211 
01212 
01213 //_____________________________________________________________________________
01214 void RooDataHist::weightError(Double_t& lo, Double_t& hi, ErrorType etype) const 
01215 { 
01216   // Return the error on current weight
01217   checkInit() ;
01218 
01219   switch (etype) {
01220 
01221   case Auto:
01222     throw string(Form("RooDataHist::weightError(%s) weight type Auto not allowed here",GetName())) ;
01223     break ;
01224 
01225   case Poisson:
01226     if (_curWgtErrLo>=0) {
01227       // Weight is preset or precalculated    
01228       lo = _curWgtErrLo ;
01229       hi = _curWgtErrHi ;
01230       return ;
01231     }
01232     
01233     // Calculate poisson errors
01234     Double_t ym,yp ;  
01235     RooHistError::instance().getPoissonInterval(Int_t(weight()+0.5),ym,yp,1) ;
01236     _curWgtErrLo = weight()-ym ;
01237     _curWgtErrHi = yp-weight() ;
01238     _errLo[_curIndex] = _curWgtErrLo ;
01239     _errHi[_curIndex] = _curWgtErrHi ;
01240     lo = _curWgtErrLo ;
01241     hi = _curWgtErrHi ;
01242     return ;
01243 
01244   case SumW2:
01245     lo = sqrt(_curSumW2) ;
01246     hi = sqrt(_curSumW2) ;
01247     return ;
01248 
01249   case None:
01250     lo = 0 ;
01251     hi = 0 ;
01252     return ;
01253   }
01254 }
01255 
01256 
01257 // wve adjust for variable bin sizes
01258 
01259 //_____________________________________________________________________________
01260 Double_t RooDataHist::interpolateDim(RooRealVar& dim, const RooAbsBinning* binning, Double_t xval, Int_t intOrder, Bool_t correctForBinSize, Bool_t cdfBoundaries) 
01261 {
01262   // Perform boundary safe 'intOrder'-th interpolation of weights in dimension 'dim'
01263   // at current value 'xval'
01264  
01265   // Fill workspace arrays spanning interpolation area
01266   Int_t fbinC = dim.getBin(*binning) ;
01267   Int_t fbinLo = fbinC-intOrder/2 - ((xval<binning->binCenter(fbinC))?1:0) ;
01268   Int_t fbinM = dim.numBins(*binning) ;
01269 
01270 
01271   Int_t i ;
01272   Double_t yarr[10] ;
01273   Double_t xarr[10] ;
01274   for (i=fbinLo ; i<=intOrder+fbinLo ; i++) {
01275     Int_t ibin ;
01276     if (i>=0 && i<fbinM) {
01277       // In range
01278       ibin = i ;
01279       dim.setBinFast(ibin,*binning) ;
01280       //cout << "INRANGE: dim.getVal(ibin=" << ibin << ") = " << dim.getVal() << endl ;
01281       xarr[i-fbinLo] = dim.getVal() ;
01282       Int_t idx = calcTreeIndex() ;      
01283       yarr[i-fbinLo] = _wgt[idx] ; 
01284       if (correctForBinSize) yarr[i-fbinLo] /=  _binv[idx] ;
01285     } else if (i>=fbinM) {
01286       // Overflow: mirror
01287       ibin = 2*fbinM-i-1 ;
01288       dim.setBinFast(ibin,*binning) ;
01289       //cout << "OVERFLOW: dim.getVal(ibin=" << ibin << ") = " << dim.getVal() << endl ;
01290       if (cdfBoundaries) {      
01291         xarr[i-fbinLo] = dim.getMax()+1e-10*(i-fbinM+1) ;
01292         yarr[i-fbinLo] = 1.0 ;
01293       } else {
01294         Int_t idx = calcTreeIndex() ;      
01295         xarr[i-fbinLo] = 2*dim.getMax()-dim.getVal() ;
01296         yarr[i-fbinLo] = _wgt[idx] ; 
01297         if (correctForBinSize) yarr[i-fbinLo] /=  _binv[idx] ;
01298       }
01299     } else {
01300       // Underflow: mirror
01301       ibin = -i - 1 ;
01302       dim.setBinFast(ibin,*binning) ;
01303       //cout << "UNDERFLOW: dim.getVal(ibin=" << ibin << ") = " << dim.getVal() << endl ;
01304       if (cdfBoundaries) {
01305         xarr[i-fbinLo] = dim.getMin()-ibin*(1e-10) ; ;
01306         yarr[i-fbinLo] = 0.0 ;
01307       } else {
01308         Int_t idx = calcTreeIndex() ;      
01309         xarr[i-fbinLo] = 2*dim.getMin()-dim.getVal() ;
01310         yarr[i-fbinLo] = _wgt[idx] ; 
01311         if (correctForBinSize) yarr[i-fbinLo] /=  _binv[idx] ;
01312       }
01313     }
01314     //cout << "ibin = " << ibin << endl ;
01315   }
01316 //   for (int k=0 ; k<=intOrder ; k++) {
01317 //     cout << "k=" << k << " x = " << xarr[k] << " y = " << yarr[k] << endl ;
01318 //   }
01319   dim.setBinFast(fbinC,*binning) ;
01320   Double_t ret = RooMath::interpolate(xarr,yarr,intOrder+1,xval) ;
01321   return ret ;
01322 }
01323 
01324 
01325 
01326 
01327 //_____________________________________________________________________________
01328 void RooDataHist::add(const RooArgSet& row, Double_t wgt, Double_t sumw2) 
01329 {
01330   // Increment the weight of the bin enclosing the coordinates given
01331   // by 'row' by the specified amount. Add the sum of weights squared
01332   // for the bin by 'sumw2' rather than wgt^2
01333 
01334   checkInit() ;
01335 
01336 //   cout << "RooDataHist::add() accepted coordinate is " << endl ;
01337 //   _vars.Print("v") ;
01338 
01339   _vars = row ;
01340   Int_t idx = calcTreeIndex() ;
01341   _wgt[idx] += wgt ;  
01342   _sumw2[idx] += (sumw2>0?sumw2:wgt*wgt) ;
01343   _errLo[idx] = -1 ;
01344   _errHi[idx] = -1 ;
01345 }
01346 
01347 
01348 
01349 //_____________________________________________________________________________
01350 void RooDataHist::set(const RooArgSet& row, Double_t wgt, Double_t wgtErrLo, Double_t wgtErrHi) 
01351 {
01352   // Increment the weight of the bin enclosing the coordinates
01353   // given by 'row' by the specified amount. Associate errors
01354   // [wgtErrLo,wgtErrHi] with the event weight on this bin.
01355 
01356   checkInit() ;
01357 
01358   _vars = row ;
01359   Int_t idx = calcTreeIndex() ;
01360   _wgt[idx] = wgt ;  
01361   _errLo[idx] = wgtErrLo ;  
01362   _errHi[idx] = wgtErrHi ;  
01363 }
01364 
01365 
01366 
01367 //_____________________________________________________________________________
01368 void RooDataHist::set(Double_t wgt, Double_t wgtErr) 
01369 {
01370   // Increment the weight of the bin enclosing the coordinates
01371   // given by 'row' by the specified amount. Associate errors
01372   // [wgtErrLo,wgtErrHi] with the event weight on this bin.
01373 
01374   checkInit() ;
01375   
01376   if (_curIndex<0) {
01377     _curIndex = calcTreeIndex() ;
01378   }
01379 
01380   _wgt[_curIndex] = wgt ;  
01381   _errLo[_curIndex] = wgtErr ;  
01382   _errHi[_curIndex] = wgtErr ;  
01383   _sumw2[_curIndex] = wgtErr*wgtErr ;
01384 }
01385 
01386 
01387 
01388 //_____________________________________________________________________________
01389 void RooDataHist::set(const RooArgSet& row, Double_t wgt, Double_t wgtErr) 
01390 {
01391   // Increment the weight of the bin enclosing the coordinates
01392   // given by 'row' by the specified amount. Associate errors
01393   // [wgtErrLo,wgtErrHi] with the event weight on this bin.
01394 
01395   checkInit() ;
01396 
01397   _vars = row ;
01398   Int_t idx = calcTreeIndex() ;
01399   _wgt[idx] = wgt ;  
01400   _errLo[idx] = wgtErr ;  
01401   _errHi[idx] = wgtErr ;  
01402   _sumw2[idx] = wgtErr*wgtErr ;
01403 }
01404 
01405 
01406 
01407 //_____________________________________________________________________________
01408 void RooDataHist::add(const RooAbsData& dset, const char* cut, Double_t wgt) 
01409 {  
01410   // Add all data points contained in 'dset' to this data set with given weight.
01411   // Optional cut string expression selects the data points to be added and can
01412   // reference any variable contained in this data set
01413 
01414   RooFormulaVar cutVar("select",cut,*dset.get()) ;
01415   add(dset,&cutVar,wgt) ;
01416 }
01417 
01418 
01419 
01420 //_____________________________________________________________________________
01421 void RooDataHist::add(const RooAbsData& dset, const RooFormulaVar* cutVar, Double_t wgt) 
01422 {
01423   // Add all data points contained in 'dset' to this data set with given weight.
01424   // Optional RooFormulaVar pointer selects the data points to be added.
01425 
01426   checkInit() ;
01427 
01428   RooFormulaVar* cloneVar = 0;
01429   RooArgSet* tmp(0) ;
01430   if (cutVar) {
01431     // Deep clone cutVar and attach clone to this dataset
01432     tmp = (RooArgSet*) RooArgSet(*cutVar).snapshot() ;
01433     if (!tmp) {
01434       coutE(DataHandling) << "RooDataHist::add(" << GetName() << ") Couldn't deep-clone cut variable, abort," << endl ;
01435       return ;
01436     }
01437 
01438     cloneVar = (RooFormulaVar*) tmp->find(cutVar->GetName()) ;
01439     cloneVar->attachDataSet(dset) ;
01440   }
01441 
01442 
01443   Int_t i ;
01444   for (i=0 ; i<dset.numEntries() ; i++) {
01445     const RooArgSet* row = dset.get(i) ;
01446     if (!cloneVar || cloneVar->getVal()) {
01447       add(*row,wgt*dset.weight()) ;
01448     }
01449   }
01450 
01451   if (cloneVar) {
01452     delete tmp ;
01453   } 
01454 }
01455 
01456 
01457 
01458 //_____________________________________________________________________________
01459 Double_t RooDataHist::sum(Bool_t correctForBinSize) const 
01460 {
01461   // Return the sum of the weights of all hist bins.
01462   //
01463   // If correctForBinSize is specified, the sum of weights
01464   // is multiplied by the N-dimensional bin volume,
01465   // making the return value the integral over the function
01466   // represented by this histogram
01467 
01468   checkInit() ;
01469     
01470   Int_t i ;
01471   Double_t total(0) ;
01472   for (i=0 ; i<_arrSize ; i++) {
01473     
01474     Double_t theBinVolume = correctForBinSize ? _binv[i] : 1.0 ;
01475     // cout << "total += " << _wgt[i] << "*" << theBinVolume << endl ;
01476     total += _wgt[i]*theBinVolume ;
01477   }
01478 
01479   return total ;
01480 }
01481 
01482 
01483 
01484 //_____________________________________________________________________________
01485 Double_t RooDataHist::sum(const RooArgSet& sumSet, const RooArgSet& sliceSet, Bool_t correctForBinSize)
01486 {
01487   // Return the sum of the weights of a multi-dimensional slice of the histogram
01488   // by summing only over the dimensions specified in sumSet.
01489   //   
01490   // The coordinates of all other dimensions are fixed to those given in sliceSet
01491   //
01492   // If correctForBinSize is specified, the sum of weights
01493   // is multiplied by the M-dimensional bin volume, (M = N(sumSet)),
01494   // making the return value the integral over the function
01495   // represented by this histogram
01496 
01497   checkInit() ;
01498 
01499   RooArgSet varSave ;
01500   varSave.addClone(_vars) ;
01501 
01502   RooArgSet* sliceOnlySet = new RooArgSet(sliceSet) ;
01503   sliceOnlySet->remove(sumSet,kTRUE,kTRUE) ;
01504 
01505   _vars = *sliceOnlySet ;
01506   calculatePartialBinVolume(*sliceOnlySet) ;
01507   delete sliceOnlySet ;
01508 
01509   TIterator* ssIter = sumSet.createIterator() ;
01510   
01511   // Calculate mask and refence plot bins for non-iterating variables
01512   RooAbsArg* arg ;
01513   Bool_t* mask = new Bool_t[_vars.getSize()] ;
01514   Int_t*  refBin = new Int_t[_vars.getSize()] ;
01515 
01516   Int_t i(0) ;
01517   _iterator->Reset() ;
01518   while((arg=(RooAbsArg*)_iterator->Next())) {
01519     if (sumSet.find(arg->GetName())) {
01520       mask[i] = kFALSE ;
01521     } else {
01522       mask[i] = kTRUE ;
01523       // coverity[FORWARD_NULL]
01524       refBin[i] = (dynamic_cast<RooAbsLValue*>(arg))->getBin() ;
01525     }
01526     i++ ;
01527   }
01528     
01529   // Loop over entire data set, skipping masked entries
01530   Double_t total(0) ;
01531   Int_t ibin ;
01532   for (ibin=0 ; ibin<_arrSize ; ibin++) {
01533 
01534     Int_t idx(0), tmp(ibin), ivar(0) ;
01535     Bool_t skip(kFALSE) ;
01536 
01537     // Check if this bin belongs in selected slice
01538     _iterator->Reset() ;
01539     // coverity[UNUSED_VALUE]
01540     while((!skip && (arg=(RooAbsArg*)_iterator->Next()))) {
01541       idx  = tmp / _idxMult[ivar] ;
01542       tmp -= idx*_idxMult[ivar] ;
01543       if (mask[ivar] && idx!=refBin[ivar]) skip=kTRUE ;
01544       ivar++ ;
01545     }
01546     
01547     if (!skip) {
01548       Double_t theBinVolume = correctForBinSize ? (*_pbinv)[ibin] : 1.0 ;
01549       //cout << "adding bin[" << ibin << "] to sum wgt = " << _wgt[ibin] << " binv = " << theBinVolume << endl ;
01550       total += _wgt[ibin]/theBinVolume ;
01551     }
01552   }
01553   delete ssIter ;
01554 
01555   delete[] mask ;
01556   delete[] refBin ;
01557 
01558   _vars = varSave ;
01559 
01560   return total ;
01561 }
01562 
01563 
01564 
01565 //_____________________________________________________________________________
01566 void RooDataHist::calculatePartialBinVolume(const RooArgSet& dimSet) const 
01567 {
01568   // Fill the transient cache with partial bin volumes with up-to-date
01569   // values for the partial volume specified by observables 'dimSet'
01570 
01571   // Allocate cache if not yet existing
01572   vector<Double_t> *pbinv = _pbinvCacheMgr.getObj(&dimSet) ;
01573   if (pbinv) {
01574     _pbinv = pbinv ;
01575     return ;
01576   }
01577 
01578   pbinv = new vector<Double_t>(_arrSize) ;
01579 
01580   // Calculate plot bins of components from master index
01581   Bool_t* selDim = new Bool_t[_vars.getSize()] ;
01582   _iterator->Reset() ;
01583   RooAbsArg* v ;
01584   Int_t i(0) ;
01585   while((v=(RooAbsArg*)_iterator->Next())) {
01586     selDim[i++] = dimSet.find(v->GetName()) ? kTRUE : kFALSE ;
01587   }
01588 
01589   // Recalculate partial bin volume cache
01590   Int_t ibin ;
01591   for (ibin=0 ; ibin<_arrSize ; ibin++) {
01592     _iterator->Reset() ;
01593     RooAbsLValue* arg ;
01594     Int_t j(0), idx(0), tmp(ibin) ;
01595     Double_t theBinVolume(1) ;
01596     while((arg=dynamic_cast<RooAbsLValue*>(_iterator->Next()))) {
01597       idx  = tmp / _idxMult[j] ;
01598       tmp -= idx*_idxMult[j++] ;
01599       if (selDim[j-1]) {
01600         RooAbsLValue* arglv = dynamic_cast<RooAbsLValue*>(arg) ;
01601         theBinVolume *= arglv->getBinWidth(idx) ;
01602       }
01603     }
01604     (*pbinv)[ibin] = theBinVolume ;
01605   }
01606 
01607   delete[] selDim ;
01608 
01609   // Put in cache (which takes ownership) 
01610   _pbinvCacheMgr.setObj(&dimSet,pbinv) ;
01611 
01612   // Publicize the array
01613   _pbinv = pbinv ;
01614 }
01615 
01616 
01617 
01618 //_____________________________________________________________________________
01619 Int_t RooDataHist::numEntries() const 
01620 {
01621   // Return the number of bins
01622 
01623   return RooAbsData::numEntries() ;
01624 }
01625 
01626 
01627 
01628 //_____________________________________________________________________________
01629 Double_t RooDataHist::sumEntries(const char* cutSpec, const char* cutRange) const
01630 {
01631   // Return the sum of weights in all entries matching cutSpec (if specified)
01632   // and in named range cutRange (if specified)
01633   // Return the
01634   checkInit() ;
01635   if (cutSpec==0 && cutRange==0) {
01636     Int_t i ;
01637     Double_t n(0) ;
01638     for (i=0 ; i<_arrSize ; i++) {
01639       if (!_binValid || _binValid[i]) {
01640         n+= _wgt[i] ;
01641       }
01642     }
01643     return n ;
01644   } else {
01645     
01646     // Setup RooFormulaVar for cutSpec if it is present
01647     RooFormula* select = 0 ;
01648     if (cutSpec) {
01649       select = new RooFormula("select",cutSpec,*get()) ;
01650     }
01651     
01652     // Otherwise sum the weights in the event
01653     Double_t sumw(0) ;
01654     Int_t i ;
01655     for (i=0 ; i<numEntries() ; i++) {
01656       get(i) ;
01657       if (select && select->eval()==0.) continue ;
01658       if (cutRange && !_vars.allInRange(cutRange)) continue ;
01659 
01660       if (!_binValid || _binValid[i]) {
01661         sumw += weight() ;
01662       }
01663     }
01664     
01665     if (select) delete select ;
01666     
01667     return sumw ;
01668   }
01669 }
01670 
01671 
01672 
01673 //_____________________________________________________________________________
01674 void RooDataHist::reset() 
01675 {
01676   // Reset all bin weights to zero
01677 
01678   // WVE DO NOT CALL RooTreeData::reset() for binned
01679   // datasets as this will delete the bin definitions
01680 
01681   Int_t i ;
01682   for (i=0 ; i<_arrSize ; i++) {
01683     _wgt[i] = 0. ;
01684     _errLo[i] = -1 ;
01685     _errHi[i] = -1 ;
01686   }
01687   _curWeight = 0 ;
01688   _curWgtErrLo = -1 ;
01689   _curWgtErrHi = -1 ;
01690   _curVolume = 1 ;
01691 
01692 }
01693 
01694 
01695 
01696 //_____________________________________________________________________________
01697 const RooArgSet* RooDataHist::get(Int_t masterIdx) const  
01698 {
01699   // Return an argset with the bin center coordinates for 
01700   // bin sequential number 'masterIdx'. For iterative use.
01701   checkInit() ;
01702   _curWeight = _wgt[masterIdx] ;
01703   _curWgtErrLo = _errLo[masterIdx] ;
01704   _curWgtErrHi = _errHi[masterIdx] ;
01705   _curSumW2 = _sumw2[masterIdx] ;
01706   _curVolume = _binv[masterIdx] ; 
01707   _curIndex  = masterIdx ;
01708   return RooAbsData::get(masterIdx) ;  
01709 }
01710 
01711 
01712 
01713 //_____________________________________________________________________________
01714 const RooArgSet* RooDataHist::get(const RooArgSet& coord) const
01715 {
01716   // Return a RooArgSet with center coordinates of the bin
01717   // enclosing the point 'coord'
01718 
01719   ((RooDataHist*)this)->_vars = coord ;
01720   return get(calcTreeIndex()) ;
01721 }
01722 
01723 
01724 
01725 //_____________________________________________________________________________
01726 Double_t RooDataHist::binVolume(const RooArgSet& coord) 
01727 {
01728   // Return the volume of the bin enclosing coordinates 'coord'
01729   checkInit() ;
01730   ((RooDataHist*)this)->_vars = coord ;
01731   return _binv[calcTreeIndex()] ;
01732 }
01733 
01734 
01735 //_____________________________________________________________________________
01736 void RooDataHist::setAllWeights(Double_t value) 
01737 {
01738   // Set all the event weight of all bins to the specified value
01739   for (Int_t i=0 ; i<_arrSize ; i++) {
01740     _wgt[i] = value ;
01741   }
01742 }
01743 
01744 
01745 
01746 //_____________________________________________________________________________
01747 TIterator* RooDataHist::sliceIterator(RooAbsArg& sliceArg, const RooArgSet& otherArgs) 
01748 {
01749   // Create an iterator over all bins in a slice defined by the subset of observables
01750   // listed in sliceArg. The position of the slice is given by otherArgs
01751 
01752   // Update to current position
01753   _vars = otherArgs ;
01754   _curIndex = calcTreeIndex() ;
01755   
01756   RooAbsArg* intArg = _vars.find(sliceArg.GetName()) ;
01757   if (!intArg) {
01758     coutE(InputArguments) << "RooDataHist::sliceIterator() variable " << sliceArg.GetName() << " is not part of this RooDataHist" << endl ;
01759     return 0 ;
01760   }
01761   return new RooDataHistSliceIter(*this,*intArg) ;
01762 }
01763 
01764 
01765 //_____________________________________________________________________________
01766 void RooDataHist::SetName(const char *name) 
01767 {
01768   // Change the name of the RooDataHist
01769 
01770   if (_dir) _dir->GetList()->Remove(this);
01771   TNamed::SetName(name) ;
01772   if (_dir) _dir->GetList()->Add(this);
01773 }
01774 
01775 
01776 //_____________________________________________________________________________
01777 void RooDataHist::SetNameTitle(const char *name, const char* title) 
01778 {
01779   // Change the title of this RooDataHist
01780 
01781   if (_dir) _dir->GetList()->Remove(this);
01782   TNamed::SetNameTitle(name,title) ;
01783   if (_dir) _dir->GetList()->Add(this);
01784 }
01785 
01786 
01787 //_____________________________________________________________________________
01788 void RooDataHist::printValue(ostream& os) const 
01789 {
01790   // Print value of the dataset, i.e. the sum of weights contained in the dataset
01791   os << numEntries() << " bins (" << sumEntries() << " weights)" ;
01792 }
01793 
01794 
01795 
01796 
01797 //_____________________________________________________________________________
01798 void RooDataHist::printArgs(ostream& os) const 
01799 {
01800   // Print argument of dataset, i.e. the observable names
01801 
01802   os << "[" ;    
01803   _iterator->Reset() ;
01804   RooAbsArg* arg ;
01805   Bool_t first(kTRUE) ;
01806   while((arg=(RooAbsArg*)_iterator->Next())) {
01807     if (first) {
01808       first=kFALSE ;
01809     } else {
01810       os << "," ;
01811     }
01812     os << arg->GetName() ;
01813   }
01814   os << "]" ;
01815 }
01816 
01817 
01818 
01819 //_____________________________________________________________________________
01820 void RooDataHist::cacheValidEntries() 
01821 {
01822   // Cache the datahist entries with bin centers that are inside/outside the
01823   // current observable definitio
01824   checkInit() ;
01825 
01826   if (!_binValid) {
01827     _binValid = new Bool_t[_arrSize] ;
01828   }
01829   TIterator* iter = _vars.createIterator() ;
01830   RooAbsArg* arg ;
01831   for (Int_t i=0 ; i<_arrSize ; i++) {
01832     get(i) ;
01833     _binValid[i] = kTRUE ;
01834     iter->Reset() ;
01835     while((arg=(RooAbsArg*)iter->Next())) {
01836       // coverity[CHECKED_RETURN]
01837       _binValid[i] &= arg->inRange(0) ;      
01838     }
01839   }
01840   delete iter ;
01841 
01842 }
01843 
01844 
01845 //_____________________________________________________________________________
01846 Bool_t RooDataHist::valid() const 
01847 {
01848   // Return true if currently loaded coordinate is considered valid within
01849   // the current range definitions of all observables
01850 
01851   // If caching is enabled, use the precached result
01852   if (_binValid) {
01853     return _binValid[_curIndex] ;
01854   }
01855 
01856   return kTRUE ;
01857 }
01858 
01859 
01860 
01861 //_____________________________________________________________________________
01862 Bool_t RooDataHist::isNonPoissonWeighted() const
01863 {
01864   // Returns true if datasets contains entries with a non-integer weight
01865 
01866   for (int i=0 ; i<numEntries() ; i++) {
01867     if (fabs(_wgt[i]-Int_t(_wgt[i]))>1e-10) return kTRUE ;
01868   }
01869   return kFALSE ;
01870 }
01871 
01872 
01873 
01874 
01875 //_____________________________________________________________________________
01876 void RooDataHist::printMultiline(ostream& os, Int_t content, Bool_t verbose, TString indent) const 
01877 {
01878   // Print the details on the dataset contents
01879 
01880   RooAbsData::printMultiline(os,content,verbose,indent) ;  
01881 
01882   os << indent << "Binned Dataset " << GetName() << " (" << GetTitle() << ")" << endl ;
01883   os << indent << "  Contains " << numEntries() << " bins with a total weight of " << sumEntries() << endl;
01884   
01885   if (!verbose) {
01886     os << indent << "  Observables " << _vars << endl ;
01887   } else {
01888     os << indent << "  Observables: " ;
01889     _vars.printStream(os,kName|kValue|kExtras|kTitle,kVerbose,indent+"  ") ;
01890   }
01891 
01892   if(verbose) {
01893     if (_cachedVars.getSize()>0) {
01894       os << indent << "  Caches " << _cachedVars << endl ;
01895     }
01896   }
01897 }
01898 
01899 
01900 
01901 //______________________________________________________________________________
01902 void RooDataHist::Streamer(TBuffer &R__b)
01903 {
01904    // Stream an object of class RooDataHist.
01905 
01906    if (R__b.IsReading()) {
01907 
01908      UInt_t R__s, R__c;
01909      Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
01910      
01911       if (R__v>2) {
01912 
01913         R__b.ReadClassBuffer(RooDataHist::Class(),this,R__v,R__s,R__c);
01914         initialize(0,kFALSE) ;
01915 
01916       } else {  
01917 
01918         // Legacy dataset conversion happens here. Legacy RooDataHist inherits from RooTreeData
01919         // which in turn inherits from RooAbsData. Manually stream RooTreeData contents on 
01920         // file here and convert it into a RooTreeDataStore which is installed in the 
01921         // new-style RooAbsData base class
01922         
01923         // --- This is the contents of the streamer code of RooTreeData version 2 ---
01924         UInt_t R__s1, R__c1;
01925         Version_t R__v1 = R__b.ReadVersion(&R__s1, &R__c1); if (R__v1) { }
01926         
01927         RooAbsData::Streamer(R__b);
01928         TTree* X_tree(0) ; R__b >> X_tree;
01929         RooArgSet X_truth ; X_truth.Streamer(R__b);
01930         TString X_blindString ; X_blindString.Streamer(R__b);
01931         R__b.CheckByteCount(R__s1, R__c1, RooTreeData::Class());
01932         // --- End of RooTreeData-v1 streamer
01933         
01934         // Construct RooTreeDataStore from X_tree and complete initialization of new-style RooAbsData
01935         _dstore = new RooTreeDataStore(X_tree,_vars) ;
01936         _dstore->SetName(GetName()) ;
01937         _dstore->SetTitle(GetTitle()) ;
01938         _dstore->checkInit() ;       
01939         
01940         RooDirItem::Streamer(R__b);
01941         R__b >> _arrSize;
01942         delete [] _wgt;
01943         _wgt = new Double_t[_arrSize];
01944         R__b.ReadFastArray(_wgt,_arrSize);
01945         delete [] _errLo;
01946         _errLo = new Double_t[_arrSize];
01947         R__b.ReadFastArray(_errLo,_arrSize);
01948         delete [] _errHi;
01949         _errHi = new Double_t[_arrSize];
01950         R__b.ReadFastArray(_errHi,_arrSize);
01951         delete [] _sumw2;
01952         _sumw2 = new Double_t[_arrSize];
01953         R__b.ReadFastArray(_sumw2,_arrSize);
01954         delete [] _binv;
01955         _binv = new Double_t[_arrSize];
01956         R__b.ReadFastArray(_binv,_arrSize);
01957         _realVars.Streamer(R__b);
01958         R__b >> _curWeight;
01959         R__b >> _curWgtErrLo;
01960         R__b >> _curWgtErrHi;
01961         R__b >> _curSumW2;
01962         R__b >> _curVolume;
01963         R__b >> _curIndex;
01964         R__b.CheckByteCount(R__s, R__c, RooDataHist::IsA());
01965 
01966       }
01967       
01968    } else {
01969 
01970       R__b.WriteClassBuffer(RooDataHist::Class(),this);
01971    }
01972 }
01973 
01974 

Generated on Tue Jul 5 15:06:25 2011 for ROOT_528-00b_version by  doxygen 1.5.1