CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_6/src/HiggsAnalysis/CombinedLimit/src/VerticalInterpHistPdf.cc

Go to the documentation of this file.
00001 #include "../interface/VerticalInterpHistPdf.h"
00002 
00003 #include <cassert>
00004 #include <memory>
00005 
00006 #include "RooFit.h"
00007 #include "Riostream.h"
00008 
00009 #include "TIterator.h"
00010 #include "RooRealVar.h"
00011 #include "RooMsgService.h"
00012 
00013 //#define TRACE_CALLS
00014 #ifdef TRACE_CALLS
00015 #include "../interface/ProfilingTools.h"
00016 #define TRACEME()   PerfCounter::add( __PRETTY_FUNCTION__ );
00017 #else
00018 #define TRACEME() 
00019 #endif
00020 
00021 ClassImp(VerticalInterpHistPdf)
00022 
00023 
00024 //_____________________________________________________________________________
00025 VerticalInterpHistPdf::VerticalInterpHistPdf() :
00026    _cacheTotal(0),
00027    _cacheSingle(0),
00028    _cacheSingleGood(0) 
00029 {
00030   // Default constructor
00031   _funcIter  = _funcList.createIterator() ;
00032   _coefIter  = _coefList.createIterator() ;
00033 }
00034 
00035 
00036 //_____________________________________________________________________________
00037 VerticalInterpHistPdf::VerticalInterpHistPdf(const char *name, const char *title, const RooRealVar &x, const RooArgList& inFuncList, const RooArgList& inCoefList, Double_t smoothRegion, Int_t smoothAlgo) :
00038   RooAbsPdf(name,title),
00039   _x("x","Independent variable",this,const_cast<RooRealVar&>(x)),
00040   _funcList("funcList","List of pdfs",this),
00041   _coefList("coefList","List of coefficients",this), // we should get shapeDirty when coefficients change
00042   _smoothRegion(smoothRegion),
00043   _smoothAlgo(smoothAlgo),
00044   _cacheTotal(0),
00045   _cacheSingle(0),
00046   _cacheSingleGood(0) 
00047 { 
00048 
00049   if (inFuncList.getSize()!=2*inCoefList.getSize()+1) {
00050     coutE(InputArguments) << "VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() 
00051                           << ") number of pdfs and coefficients inconsistent, must have Nfunc=1+2*Ncoef" << endl ;
00052     assert(0);
00053   }
00054 
00055   TIterator* funcIter = inFuncList.createIterator() ;
00056   RooAbsArg* func;
00057   while((func = (RooAbsArg*)funcIter->Next())) {
00058     RooAbsPdf *pdf = dynamic_cast<RooAbsPdf*>(func);
00059     if (!pdf) {
00060       coutE(InputArguments) << "ERROR: VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() << ") function  " << func->GetName() << " is not of type RooAbsPdf" << endl;
00061       assert(0);
00062     }
00063     RooArgSet *params = pdf->getParameters(RooArgSet(x));
00064     if (params->getSize() > 0) {
00065       coutE(InputArguments) << "ERROR: VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() << ") pdf  " << func->GetName() << " has some parameters." << endl;
00066       assert(0);
00067     }
00068     delete params;
00069     _funcList.add(*func) ;
00070   }
00071   delete funcIter;
00072 
00073   TIterator* coefIter = inCoefList.createIterator() ;
00074   RooAbsArg* coef;
00075   while((coef = (RooAbsArg*)coefIter->Next())) {
00076     if (!dynamic_cast<RooAbsReal*>(coef)) {
00077       coutE(InputArguments) << "ERROR: VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() << ") coefficient " << coef->GetName() << " is not of type RooAbsReal" << endl;
00078       assert(0);
00079     }
00080     _coefList.add(*coef) ;    
00081   }
00082   delete coefIter;
00083 
00084   _funcIter  = _funcList.createIterator() ;
00085   _coefIter = _coefList.createIterator() ;
00086 
00087 }
00088 
00089 
00090 
00091 
00092 //_____________________________________________________________________________
00093 VerticalInterpHistPdf::VerticalInterpHistPdf(const VerticalInterpHistPdf& other, const char* name) :
00094   RooAbsPdf(other,name),
00095   _x("x",this,other._x),
00096   _funcList("funcList",this,other._funcList),
00097   _coefList("coefList",this,other._coefList),
00098   _smoothRegion(other._smoothRegion),
00099   _smoothAlgo(other._smoothAlgo),
00100   _cacheTotal(0),
00101   _cacheSingle(0),
00102   _cacheSingleGood(0) 
00103 {
00104   // Copy constructor
00105 
00106   _funcIter  = _funcList.createIterator() ;
00107   _coefIter = _coefList.createIterator() ;
00108 }
00109 
00110 
00111 
00112 //_____________________________________________________________________________
00113 VerticalInterpHistPdf::~VerticalInterpHistPdf()
00114 {
00115   // Destructor
00116   delete _funcIter ;
00117   delete _coefIter ;
00118   if (_cacheTotal) {
00119       delete _cacheTotal;
00120       for (int i = 0; i < _funcList.getSize(); ++i) delete _cacheSingle[i]; 
00121       delete [] _cacheSingle;
00122       delete [] _cacheSingleGood;
00123   }
00124 }
00125 
00126 //_____________________________________________________________________________
00127 Double_t VerticalInterpHistPdf::evaluate() const 
00128 {
00129   if (_cacheTotal == 0) setupCaches();
00130 #if 0
00131   printf("Evaluate called for x %f, cache status %d: ", _x.arg().getVal(), _sentry.good());
00132   int ndim = _coefList.getSize();
00133   for (int i = 0; i < ndim; ++i) { 
00134     RooAbsReal *rar = dynamic_cast<RooAbsReal *>(_coefList.at(i));
00135     printf("%s = %+6.4f  ", rar->GetName(), rar->getVal());
00136   }
00137   std::cout << std::endl;
00138 #endif
00139 
00140   if (!_sentry.good()) syncTotal();
00141   int nbin = _cacheTotal->GetNbinsX(), ibin = _cacheTotal->FindBin(_x);
00142   if (ibin < 1) ibin = 1;
00143   else if (ibin > nbin) ibin = nbin;
00144   return _cacheTotal->GetBinContent(ibin);
00145 }
00146 
00147 
00148 void VerticalInterpHistPdf::syncComponent(int i) const {
00149     RooAbsPdf *pdfi = dynamic_cast<RooAbsPdf *>(_funcList.at(i));
00150     if (_cacheSingle[i] != 0) delete _cacheSingle[i];
00151     _cacheSingle[i] = pdfi->createHistogram("",dynamic_cast<const RooRealVar &>(_x.arg())); 
00152     _cacheSingle[i]->SetDirectory(0);
00153     if (_cacheSingle[i]->Integral("width")) { _cacheSingle[i]->Scale(1.0/_cacheSingle[i]->Integral("width")); }
00154     if (i > 0) {
00155         for (int b = 1, nb = _cacheSingle[i]->GetNbinsX(); b <= nb; ++b) {
00156             double y  = _cacheSingle[i]->GetBinContent(b);
00157             double y0 = _cacheSingle[0]->GetBinContent(b);
00158             if (_smoothAlgo < 0) {
00159                 if (y > 0 && y0 > 0) {
00160                     double logk = log(y/y0);
00161                     // odd numbers correspond to up variations, even numbers to down variations,
00162                     // and down variations need -log(kappa) instead of log(kappa)
00163                     _cacheSingle[i]->SetBinContent(b, logk);
00164                 } else {
00165                     _cacheSingle[i]->SetBinContent(b, 0);
00166                 }
00167             } else {
00168                 _cacheSingle[i]->SetBinContent(b, y - y0);
00169             }
00170         }
00171     }
00172     _cacheSingleGood[i] = true;
00173 }
00174 
00175 void VerticalInterpHistPdf::syncTotal() const {
00176     int ndim = _coefList.getSize();
00177     for (int i = 0; i < 2*ndim+1; ++i) {
00178         if (!_cacheSingleGood[i]) syncComponent(i);
00179     }
00180     for (int b = 1, nb = _cacheTotal->GetNbinsX(); b <= nb; ++b) {
00181         double val = _cacheSingle[0]->GetBinContent(b);
00182         _coefIter->Reset();
00183         for (int i = 0; i < ndim; ++i) {
00184             double dhi = _cacheSingle[2*i+1]->GetBinContent(b);
00185             double dlo = _cacheSingle[2*i+2]->GetBinContent(b);
00186             double x = (dynamic_cast<RooAbsReal *>(_coefIter->Next()))->getVal();
00187             double alpha = x * 0.5 * ((dhi-dlo) + (dhi+dlo)*smoothStepFunc(x));
00188             // alpha(0) = 0
00189             // alpha(+1) = dhi 
00190             // alpha(-1) = dlo
00191             // alpha(x >= +1) = |x|*dhi
00192             // alpha(x <= -1) = |x|*dlo
00193             // alpha is continuous and has continuous first and second derivative, as smoothStepFunc has them
00194             if (_smoothAlgo < 0) {
00195                 val *= exp(alpha);
00196             } else {
00197                 val += alpha; 
00198             }
00199         }    
00200         if (val <= 0) val = 1e-9;
00201         _cacheTotal->SetBinContent(b, val);
00202     }
00203     double norm = _cacheTotal->Integral("width");
00204     if (norm > 0) _cacheTotal->Scale(1.0/norm);
00205     _sentry.reset();
00206 }
00207 
00208 void VerticalInterpHistPdf::setupCaches() const {
00209     int ndim = _coefList.getSize();
00210     _cacheTotal     = (dynamic_cast<const RooRealVar &>(_x.arg())).createHistogram("total"); 
00211     _cacheTotal->SetDirectory(0);
00212     _cacheSingle     = new TH1*[2*ndim+1]; assert(_cacheSingle);
00213     _cacheSingleGood = new int[2*ndim+1];  assert(_cacheSingleGood);
00214     for (int i = 0; i < 2*ndim+1; ++i) { 
00215         _cacheSingle[i] = 0; 
00216         _cacheSingleGood[i] = 0; 
00217         syncComponent(i);  
00218     } 
00219     if (_sentry.empty()) _sentry.addVars(_coefList); 
00220     syncTotal();
00221 }
00222 
00223 //=============================================================================================
00224 ClassImp(FastVerticalInterpHistPdfBase)
00225 ClassImp(FastVerticalInterpHistPdf)
00226 ClassImp(FastVerticalInterpHistPdf2D)
00227 
00228 
00229 //_____________________________________________________________________________
00230 FastVerticalInterpHistPdfBase::FastVerticalInterpHistPdfBase() 
00231 {
00232   // Default constructor
00233   _funcIter  = _funcList.createIterator() ;
00234   _coefIter  = _coefList.createIterator() ;
00235 }
00236 
00237 
00238 //_____________________________________________________________________________
00239 FastVerticalInterpHistPdfBase::FastVerticalInterpHistPdfBase(const char *name, const char *title, const RooArgSet &obs, const RooArgList& inFuncList, const RooArgList& inCoefList, Double_t smoothRegion, Int_t smoothAlgo) :
00240   RooAbsPdf(name,title),
00241   _funcList("funcList","List of pdfs",this),
00242   _coefList("coefList","List of coefficients",this), // we should get shapeDirty when coefficients change
00243   _smoothRegion(smoothRegion),
00244   _smoothAlgo(smoothAlgo),
00245   _init(false),
00246   _morphs(), _morphParams()
00247 { 
00248   TRACEME()
00249 
00250   if (inFuncList.getSize()!=2*inCoefList.getSize()+1) {
00251     coutE(InputArguments) << "VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() 
00252                           << ") number of pdfs and coefficients inconsistent, must have Nfunc=1+2*Ncoef" << endl ;
00253     assert(0);
00254   }
00255 
00256   TIterator* funcIter = inFuncList.createIterator() ;
00257   RooAbsArg* func;
00258   while((func = (RooAbsArg*)funcIter->Next())) {
00259     RooAbsPdf *pdf = dynamic_cast<RooAbsPdf*>(func);
00260     if (!pdf) {
00261       coutE(InputArguments) << "ERROR: VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() << ") function  " << func->GetName() << " is not of type RooAbsPdf" << endl;
00262       assert(0);
00263     }
00264     RooArgSet *params = pdf->getParameters(obs);
00265     if (params->getSize() > 0) {
00266       coutE(InputArguments) << "ERROR: VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() << ") pdf  " << func->GetName() << " has some parameters." << endl;
00267       assert(0);
00268     }
00269     delete params;
00270     _funcList.add(*func) ;
00271   }
00272   delete funcIter;
00273 
00274   TIterator* coefIter = inCoefList.createIterator() ;
00275   RooAbsArg* coef;
00276   while((coef = (RooAbsArg*)coefIter->Next())) {
00277     if (!dynamic_cast<RooAbsReal*>(coef)) {
00278       coutE(InputArguments) << "ERROR: VerticalInterpHistPdf::VerticalInterpHistPdf(" << GetName() << ") coefficient " << coef->GetName() << " is not of type RooAbsReal" << endl;
00279       assert(0);
00280     }
00281     _coefList.add(*coef) ;    
00282   }
00283   delete coefIter;
00284 
00285   _funcIter  = _funcList.createIterator() ;
00286   _coefIter = _coefList.createIterator() ;
00287 
00288 }
00289 
00290 //_____________________________________________________________________________
00291 FastVerticalInterpHistPdfBase::FastVerticalInterpHistPdfBase(const FastVerticalInterpHistPdfBase& other, const char* name) :
00292   RooAbsPdf(other,name),
00293   _funcList("funcList",this,other._funcList),
00294   _coefList("coefList",this,other._coefList),
00295   _smoothRegion(other._smoothRegion),
00296   _smoothAlgo(other._smoothAlgo),
00297   _init(false),
00298 #ifdef FastVerticalInterpHistPdf_CopyConstructor
00299   _morphs(other._morphs), _morphParams(other._morphParams)
00300 #else
00301   _morphs(), _morphParams()
00302 #endif
00303 {
00304   // Copy constructor
00305   _funcIter  = _funcList.createIterator() ;
00306   _coefIter = _coefList.createIterator() ;
00307 #ifdef FastVerticalInterpHistPdf_CopyConstructor
00308   _sentry.addVars(_coefList);
00309   _sentry.setValueDirty(); 
00310 #endif
00311 }
00312 
00313 
00314 //_____________________________________________________________________________
00315 FastVerticalInterpHistPdfBase::~FastVerticalInterpHistPdfBase()
00316 {
00317   // Destructor
00318   delete _funcIter ;
00319   delete _coefIter ;
00320 }
00321 
00322 
00323 //_____________________________________________________________________________
00324 Double_t FastVerticalInterpHistPdf::evaluate() const 
00325 {
00326   TRACEME()
00327   if (_cache.size() == 0) setupCaches();
00328   if (!_sentry.good() || !_init) syncTotal();
00329   //return std::max<double>(1e-9, _cache.GetAt(_x));
00330   return _cache.GetAt(_x);
00331 }
00332 
00333 //_____________________________________________________________________________
00334 Double_t FastVerticalInterpHistPdf2D::evaluate() const 
00335 {
00336   TRACEME()
00337   if (_cache.size() == 0) setupCaches();
00338 
00339   if (!_sentry.good() || !_init) syncTotal();
00340   //return std::max<double>(1e-9, _cache.GetAt(_x, _y));
00341   return _cache.GetAt(_x, _y);
00342 }
00343 
00344 
00345 
00346 void FastVerticalInterpHistPdf::syncNominal() const {
00347     TRACEME()
00348     RooAbsPdf *pdf = dynamic_cast<RooAbsPdf *>(_funcList.at(0));
00349     const RooRealVar &x = dynamic_cast<const RooRealVar &>(_x.arg());
00350     std::auto_ptr<TH1> hist(pdf->createHistogram("",x));
00351     hist->SetDirectory(0); 
00352     _cacheNominal = FastHisto(*hist);
00353     _cacheNominal.Normalize();
00354     if (_smoothAlgo < 0) {
00355         _cacheNominalLog = _cacheNominal;
00356         _cacheNominalLog.Log();
00357     }
00358 }
00359 
00360 void FastVerticalInterpHistPdf2D::syncNominal() const {
00361     TRACEME()
00362     RooAbsPdf *pdf = dynamic_cast<RooAbsPdf *>(_funcList.at(0));
00363     const RooRealVar &x = dynamic_cast<const RooRealVar &>(_x.arg());
00364     const RooRealVar &y = dynamic_cast<const RooRealVar &>(_y.arg());
00365     const RooCmdArg &cond = _conditional ? RooFit::ConditionalObservables(RooArgSet(x)) : RooCmdArg::none();
00366     std::auto_ptr<TH1> hist(pdf->createHistogram("", x, RooFit::YVar(y), cond));
00367     hist->SetDirectory(0); 
00368     _cacheNominal = FastHisto2D(dynamic_cast<TH2F&>(*hist), _conditional);
00369     if (_conditional) _cacheNominal.NormalizeXSlices(); 
00370     else              _cacheNominal.Normalize(); 
00371 
00372     if (_smoothAlgo < 0) {
00373         _cacheNominalLog = _cacheNominal;
00374         _cacheNominalLog.Log();
00375     }
00376 }
00377 
00378 
00379 
00380 void FastVerticalInterpHistPdfBase::syncMorph(Morph &out, const FastTemplate &nominal, FastTemplate &lo, FastTemplate &hi) const {
00381     if (_smoothAlgo < 0)  {
00382         hi.LogRatio(nominal); lo.LogRatio(nominal);
00383         //printf("Log-ratios for dimension %d: \n", dim);  hi.Dump(); lo.Dump();
00384     } else {
00385         hi.Subtract(nominal); lo.Subtract(nominal);
00386         //printf("Differences for dimension %d: \n", dim);  hi.Dump(); lo.Dump();
00387     }
00388     FastTemplate::SumDiff(hi, lo, out.sum, out.diff);
00389     //printf("Sum and diff for dimension %d: \n", dim);  out.sum.Dump(); out.diff.Dump();
00390 }
00391 
00392 void FastVerticalInterpHistPdf::syncComponents(int dim) const {
00393     TRACEME()
00394     RooAbsPdf *pdfHi = dynamic_cast<RooAbsPdf *>(_funcList.at(2*dim+1));
00395     RooAbsPdf *pdfLo = dynamic_cast<RooAbsPdf *>(_funcList.at(2*dim+2));
00396     const RooRealVar &x = dynamic_cast<const RooRealVar &>(_x.arg());
00397     std::auto_ptr<TH1> histHi(pdfHi->createHistogram("",x)); histHi->SetDirectory(0); 
00398     std::auto_ptr<TH1> histLo(pdfLo->createHistogram("",x)); histLo->SetDirectory(0);
00399 
00400     FastHisto hi(*histHi), lo(*histLo); 
00401     //printf("Un-normalized templates for dimension %d: \n", dim);  hi.Dump(); lo.Dump();
00402     hi.Normalize(); lo.Normalize();
00403     //printf("Normalized templates for dimension %d: \n", dim);  hi.Dump(); lo.Dump();
00404     syncMorph(_morphs[dim], _cacheNominal, lo, hi);
00405 }
00406 void FastVerticalInterpHistPdf2D::syncComponents(int dim) const {
00407     RooAbsPdf *pdfHi = dynamic_cast<RooAbsPdf *>(_funcList.at(2*dim+1));
00408     RooAbsPdf *pdfLo = dynamic_cast<RooAbsPdf *>(_funcList.at(2*dim+2));
00409     const RooRealVar &x = dynamic_cast<const RooRealVar &>(_x.arg());
00410     const RooRealVar &y = dynamic_cast<const RooRealVar &>(_y.arg());
00411     const RooCmdArg &cond = _conditional ? RooFit::ConditionalObservables(RooArgSet(x)) : RooCmdArg::none();
00412     std::auto_ptr<TH1> histHi(pdfHi->createHistogram("", x, RooFit::YVar(y), cond)); histHi->SetDirectory(0); 
00413     std::auto_ptr<TH1> histLo(pdfLo->createHistogram("", x, RooFit::YVar(y), cond)); histLo->SetDirectory(0);
00414 
00415     FastHisto2D hi(dynamic_cast<TH2&>(*histHi), _conditional), lo(dynamic_cast<TH2&>(*histLo), _conditional); 
00416     //printf("Un-normalized templates for dimension %d: \n", dim);  hi.Dump(); lo.Dump();
00417     if (_conditional) {
00418         hi.NormalizeXSlices(); lo.NormalizeXSlices();
00419     } else {
00420         hi.Normalize(); lo.Normalize();
00421     }
00422     //printf("Normalized templates for dimension %d: \n", dim);  hi.Dump(); lo.Dump();
00423     syncMorph(_morphs[dim], _cacheNominal, lo, hi);
00424 }
00425 
00426 
00427 void FastVerticalInterpHistPdfBase::syncTotal(FastTemplate &cache, const FastTemplate &cacheNominal, const FastTemplate &cacheNominalLog) const {
00428     TRACEME()
00429     /* === how the algorithm works, in theory ===
00430      * let  dhi = h_hi - h_nominal
00431      *      dlo = h_lo - h_nominal
00432      * and x be the morphing parameter
00433      * we define alpha = x * 0.5 * ((dhi-dlo) + (dhi+dlo)*smoothStepFunc(x));
00434      * which satisfies:
00435      *     alpha(0) = 0
00436      *     alpha(+1) = dhi 
00437      *     alpha(-1) = dlo
00438      *     alpha(x >= +1) = |x|*dhi
00439      *     alpha(x <= -1) = |x|*dlo
00440      *     alpha is continuous and has continuous first and second derivative, as smoothStepFunc has them
00441      * === and in practice ===
00442      * we already have computed the histogram for diff=(dhi-dlo) and sum=(dhi+dlo)
00443      * so we just do template += (0.5 * x) * (diff + smoothStepFunc(x) * sum)
00444      * ========================================== */
00445 
00446     // start from nominal
00447     cache.CopyValues(_smoothAlgo < 0 ? cacheNominalLog : cacheNominal);
00448     //printf("Cache initialized to nominal template: \n");  cacheNominal.Dump();
00449 
00450     // apply all morphs one by one
00451     for (int i = 0, ndim = _coefList.getSize(); i < ndim; ++i) {
00452         double x = _morphParams[i]->getVal();
00453         double a = 0.5*x, b = smoothStepFunc(x);
00454         cache.Meld(_morphs[i].diff, _morphs[i].sum, a, b);    
00455         //printf("Merged transformation for dimension %d, x = %+5.3f, step = %.3f: \n", i, x, b);  cache.Dump();
00456     }
00457 
00458     // if necessary go back to linear scale
00459     if (_smoothAlgo < 0) {
00460         cache.Exp();
00461         //printf("Done exponential tranformation\n");  cache.Dump();
00462     } else {
00463         cache.CropUnderflows();
00464     }
00465     
00466     // mark as done
00467     _sentry.reset();
00468     _init = true;
00469 }
00470 
00471 void FastVerticalInterpHistPdf::syncTotal() const {
00472     FastVerticalInterpHistPdfBase::syncTotal(_cache, _cacheNominal, _cacheNominalLog);
00473 
00474     // normalize the result
00475     _cache.Normalize(); 
00476     //printf("Normalized result\n");  _cache.Dump();
00477 }
00478 
00479 void FastVerticalInterpHistPdf2D::syncTotal() const {
00480     FastVerticalInterpHistPdfBase::syncTotal(_cache, _cacheNominal, _cacheNominalLog);
00481     // normalize the result
00482     if (_conditional) _cache.NormalizeXSlices(); 
00483     else              _cache.Normalize(); 
00484     //printf("Normalized result\n");  _cache.Dump();
00485 }
00486 
00487 
00488 
00489 void FastVerticalInterpHistPdf::setupCaches() const {
00490     TRACEME()
00491     int ndim = _coefList.getSize();
00492 
00493     _morphs.resize(ndim);
00494     _morphParams.resize(ndim);
00495     syncNominal();
00496     //printf("Nominal template has been set up: \n");  _cacheNominal.Dump();
00497     _coefIter->Reset();
00498     for (int i = 0; i < ndim; ++i) {
00499         _morphParams[i] = dynamic_cast<RooAbsReal *>(_coefIter->Next());
00500         _morphs[i].sum.Resize(_cacheNominal.size());
00501         _morphs[i].diff.Resize(_cacheNominal.size());
00502         syncComponents(i);
00503     } 
00504     _cache = FastHisto(_cacheNominal);
00505 
00506     if (_sentry.empty()) _sentry.addVars(_coefList); 
00507     syncTotal();
00508 }
00509 
00510 void FastVerticalInterpHistPdf2D::setupCaches() const {
00511     TRACEME()
00512     int ndim = _coefList.getSize();
00513 
00514     _morphs.resize(ndim);
00515     _morphParams.resize(ndim);
00516     syncNominal();
00517     //printf("Nominal template has been set up: \n");  _cacheNominal.Dump();
00518     _coefIter->Reset();
00519     for (int i = 0; i < ndim; ++i) {
00520         _morphParams[i] = dynamic_cast<RooAbsReal *>(_coefIter->Next());
00521         _morphs[i].sum.Resize(_cacheNominal.size());
00522         _morphs[i].diff.Resize(_cacheNominal.size());
00523         syncComponents(i);
00524     } 
00525     _cache = FastHisto2D(_cacheNominal);
00526 
00527     if (_sentry.empty()) _sentry.addVars(_coefList); 
00528     syncTotal();
00529 }
00530 
00531