42 : _cacheMgr(this, 10, true, true), _curNormSet(0), _M(0), _MSqr(0), _setting(
RooMomentMorphFuncND::
Linear), _useHorizMorph(true)
51 :
RooAbsReal(
name, title), _cacheMgr(this, 10, true, true), _parList(
"parList",
"List of morph parameters", this),
52 _obsList(
"obsList",
"List of observables", this), _referenceGrid(referenceGrid),
53 _pdfList(
"pdfList",
"List of pdfs", this), _setting(setting), _useHorizMorph(true)
73 :
RooAbsReal(
name, title), _cacheMgr(this, 10, true, true), _parList(
"parList",
"List of morph parameters", this),
74 _obsList(
"obsList",
"List of observables", this), _pdfList(
"pdfList",
"List of pdfs", this), _setting(setting),
81 for (
int i = 0; i < mrefpoints.
GetNrows(); ++i) {
83 if (mrefpoints[i] == grid.
array()[j]) {
110 :
RooAbsReal(
name, title), _cacheMgr(this, 10, true, true), _parList(
"parList",
"List of morph parameters", this),
111 _obsList(
"obsList",
"List of observables", this), _pdfList(
"pdfList",
"List of pdfs", this), _setting(setting),
117 for (
auto *mref : mrefList) {
120 <<
" is not of type RooAbsReal" << endl;
121 throw string(
"RooMomentMorphFuncND::ctor() ERROR mref is not of type RooAbsReal");
125 <<
" is not a constant, taking a snapshot of its value" << endl;
133 for (i = 0; i < mrefpoints.
GetNrows(); ++i) {
135 if (mrefpoints[i] == grid.
array()[j]) {
160 :
RooAbsReal(other,
name), _cacheMgr(other._cacheMgr, this), _curNormSet(0),
161 _parList(
"parList", this, other._parList), _obsList(
"obsList", this, other._obsList),
162 _referenceGrid(other._referenceGrid), _pdfList(
"pdfList", this, other._pdfList), _M(0), _MSqr(0),
163 _setting(other._setting), _useHorizMorph(other._useHorizMorph)
190 for (
auto *par : parList) {
193 <<
" is not of type RooAbsReal" << endl;
194 throw string(
"RooMomentMorphFuncND::initializeParameters() ERROR parameter is not of type RooAbsReal");
203 for (
auto *var : obsList) {
206 <<
" is not of type RooAbsReal" << endl;
207 throw string(
"RooMomentMorphFuncND::initializeObservables() ERROR variable is not of type RooAbsReal");
218 vector<Digits<T>> vd;
220 for (
typename vector<vector<T>>::const_iterator it = in.begin(); it != in.end(); ++it) {
221 Digits<T> d = {(*it).begin(), (*it).end(), (*it).begin()};
227 for (
typename vector<
Digits<T>>::const_iterator it = vd.
begin(); it != vd.end(); ++it) {
228 result.push_back(*(it->me));
234 if (it->me == it->end) {
235 if (it + 1 == vd.end()) {
264 <<
": " << nPar <<
" !=" << nDim << endl;
270 <<
": " << nPdf <<
" !=" << nRef << endl;
280 vector<vector<double>> dm(nPdf);
281 for (
int k = 0; k < nPdf; ++k) {
283 for (
int idim = 0; idim < nPar; idim++) {
285 dm2.push_back(delta);
290 vector<vector<int>> powers;
291 for (
int idim = 0; idim < nPar; idim++) {
296 powers.push_back(xtmp);
299 vector<vector<int>>
output;
301 int nCombs =
output.size();
303 for (
int k = 0; k < nPdf; ++k) {
305 for (
int i = 0; i < nCombs; i++) {
307 for (
int ix = 0; ix < nPar; ix++) {
308 double delta = dm[k][ix];
327 : _grid(other._grid), _pdfList(other._pdfList), _pdfMap(other._pdfMap), _nref(other._nref)
339 vector<int> thisBoundaries;
340 vector<double> thisBoundaryCoordinates;
341 thisBoundaries.push_back(bin_x);
342 thisBoundaryCoordinates.push_back(_grid[0]->array()[bin_x]);
345 _nref.push_back(thisBoundaryCoordinates);
351 vector<int> thisBoundaries;
352 vector<double> thisBoundaryCoordinates;
353 thisBoundaries.push_back(bin_x);
354 thisBoundaryCoordinates.push_back(_grid[0]->array()[bin_x]);
355 thisBoundaries.push_back(bin_y);
356 thisBoundaryCoordinates.push_back(_grid[1]->array()[bin_y]);
359 _nref.push_back(thisBoundaryCoordinates);
365 vector<int> thisBoundaries;
366 vector<double> thisBoundaryCoordinates;
367 thisBoundaries.push_back(bin_x);
368 thisBoundaryCoordinates.push_back(_grid[0]->array()[bin_x]);
369 thisBoundaries.push_back(bin_y);
370 thisBoundaryCoordinates.push_back(_grid[1]->array()[bin_y]);
371 thisBoundaries.push_back(bin_z);
372 thisBoundaryCoordinates.push_back(_grid[2]->array()[bin_z]);
375 _nref.push_back(thisBoundaryCoordinates);
381 vector<double> thisBoundaryCoordinates;
382 int nBins = bins.size();
383 for (
int i = 0; i < nBins; i++) {
384 thisBoundaryCoordinates.push_back(_grid[i]->array()[bins[i]]);
388 _nref.push_back(thisBoundaryCoordinates);
403 vector<RooAbsReal *> meanrv(nPdf * nObs,
null);
404 vector<RooAbsReal *> sigmarv(nPdf * nObs,
null);
405 vector<RooAbsReal *> myrms(nObs,
null);
406 vector<RooAbsReal *> mypos(nObs,
null);
407 vector<RooAbsReal *> slope(nPdf * nObs,
null);
408 vector<RooAbsReal *> offsets(nPdf * nObs,
null);
409 vector<RooAbsReal *> transVar(nPdf * nObs,
null);
410 vector<RooAbsReal *> transPdf(nPdf,
null);
420 for (
int i = 0; i < 3 * nPdf; ++i) {
421 string fracName =
Form(
"frac_%d", i);
422 double initval = 0.0;
428 else if (i < 2 * nPdf)
436 string sumfuncName =
Form(
"%s_sumfunc",
GetName());
441 for (
int i = 0; i < nPdf; ++i) {
442 for (
int j = 0; j < nObs; ++j) {
450 sigmarv[
sij(i, j)] = mom;
451 meanrv[
sij(i, j)] = mom->
mean();
453 ownedComps.
add(*sigmarv[
sij(i, j)]);
458 for (
int j = 0; j < nObs; ++j) {
461 for (
int i = 0; i < nPdf; ++i) {
462 meanList.
add(*meanrv[
sij(i, j)]);
463 rmsList.
add(*sigmarv[
sij(i, j)]);
465 string myrmsName =
Form(
"%s_rms_%d",
GetName(), j);
466 string myposName =
Form(
"%s_pos_%d",
GetName(), j);
467 mypos[j] =
new RooAddition(myposName.c_str(), myposName.c_str(), meanList, coefList2);
468 myrms[j] =
new RooAddition(myrmsName.c_str(), myrmsName.c_str(), rmsList, coefList3);
476 for (
auto const *pdf : static_range_cast<RooAbsReal *>(
_pdfList)) {
478 string pdfName =
Form(
"pdf_%d", i);
482 for (
auto *var : static_range_cast<RooRealVar *>(obsList)) {
484 string slopeName =
Form(
"%s_slope_%d_%d",
GetName(), i, j);
485 string offsetName =
Form(
"%s_offset_%d_%d",
GetName(), i, j);
494 string transVarName =
Form(
"%s_transVar_%d_%d",
GetName(), i, j);
495 transVar[
sij(i, j)] =
new RooLinearVar(transVarName.c_str(), transVarName.c_str(), *var, *slope[
sij(i, j)],
496 *offsets[
sij(i, j)]);
502 ownedComps.
add(*transVar[
sij(i, j)]);
507 transPdfList.
add(*transPdf[i]);
508 ownedComps.
add(*transPdf[i]);
513 theSumFunc =
new RooRealSumFunc(sumfuncName.c_str(), sumfuncName.c_str(), transPdfList, coefList);
524 string trackerName =
Form(
"%s_frac_tracker",
GetName());
528 cache =
new CacheElem(*theSumFunc, *tracker, fracl);
594template <
typename Iterator>
597 if ((
first == last) || (
first == k) || (last == k)) {
600 Iterator itr1 =
first;
601 Iterator itr2 = last;
610 while (
first != itr1) {
611 if (*--itr1 < *itr2) {
613 while (!(*itr1 < *j)) ++j;
618 rotate(itr1, j, last);
623 rotate(k, itr2, last);
627 rotate(
first, k, last);
637 double fracLinear(1.);
638 double fracNonLinear(1.);
643 for (
int idim = 0; idim < nPar; idim++) {
645 dm2.push_back(delta);
648 vector<vector<int>> powers;
649 for (
int idim = 0; idim < nPar; idim++) {
654 powers.push_back(xtmp);
657 vector<vector<int>>
output;
659 int nCombs =
output.size();
661 vector<double> deltavec(nPdf, 1.0);
664 for (
int i = 0; i < nCombs; i++) {
666 for (
int ix = 0; ix < nPar; ix++) {
667 double delta = dm2[ix];
670 deltavec[nperm] = tmpDm;
674 double sumposfrac = 0.0;
675 for (
int i = 0; i < nPdf; ++i) {
678 for (
int j = 0; j < nPdf; ++j) {
679 ffrac += (*self.
_M)(j, i) * deltavec[j] * fracNonLinear;
692 ((
RooRealVar *)frac(nPdf + i))->setVal(ffrac);
693 ((
RooRealVar *)frac(2 * nPdf + i))->setVal(ffrac);
696 cout <<
"NonLinear fraction " << ffrac << endl;
698 frac(nPdf + i)->Print();
699 frac(2 * nPdf + i)->Print();
704 for (
int i = 0; i < nPdf; ++i) {
718 for (
int i = 0; i < nPdf; ++i) {
721 ((
RooRealVar *)frac(nPdf + i))->setVal(initval);
722 ((
RooRealVar *)frac(2 * nPdf + i))->setVal(initval);
727 for (
int j = 0; j < nPar; j++) {
729 mtmp.push_back(
m->getVal());
735 vector<double> deltavec(depth, 1.0);
740 for (
int ix = 0; ix < nPar; ix++) {
744 for (
int iperm = 1; iperm <= nPar; ++iperm) {
746 double dtmp = mtmp[xtmp[0]] - self.
_squareVec[0][xtmp[0]];
747 for (
int itmp = 1; itmp < iperm; ++itmp) {
748 dtmp *= mtmp[xtmp[itmp]] - self.
_squareVec[0][xtmp[itmp]];
750 deltavec[nperm + 1] = dtmp;
755 double origFrac1(0.), origFrac2(0.);
756 for (
int i = 0; i < depth; ++i) {
758 for (
int j = 0; j < depth; ++j) {
759 ffrac += (*self.
_MSqr)(j, i) * deltavec[j] * fracLinear;
775 cout <<
"Linear fraction " << ffrac << endl;
791 bool isEnclosed =
true;
792 for (
int i = 0; i < nPar; i++) {
803 vector<vector<double>> boundaries(nPar);
804 for (
int idim = 0; idim < nPar; idim++) {
808 boundaries[idim].push_back(lo);
809 boundaries[idim].push_back(
hi);
812 vector<vector<double>>
output;
816 for (
int isq = 0; isq < depth; isq++) {
817 for (
int iref = 0; iref < nRef; iref++) {
840 for (
int ix = 0; ix < nPar; ix++) {
844 for (
int k = 0; k < depth; ++k) {
850 for (
int iperm = 1; iperm <= nPar; ++iperm) {
852 double dtmp =
_squareVec[k][xtmp[0]] - squareBase[xtmp[0]];
853 for (
int itmp = 1; itmp < iperm; ++itmp) {
854 dtmp *=
_squareVec[k][xtmp[itmp]] - squareBase[xtmp[itmp]];
856 M(k, nperm + 1) = dtmp;
876 cout <<
"Currently BinIntegrator only knows how to deal with 1-d " << endl;
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
TMatrixT< Double_t > TMatrixD
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
void setLocalNoDirtyInhibit(bool flag) const
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
void addServerList(RooAbsCollection &serverList, bool valueProp=true, bool shapeProp=false)
Register a list of RooAbsArg as servers to us by calling addServer() for each arg in the list.
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
bool setRealValue(const char *name, double newVal=0, bool verbose=false)
Set value of a RooAbsRealLValye stored in set with given name to newVal No error messages are printed...
RooAbsArg * first() const
RooAbsMoment represents the first, second, or third order derivative of any RooAbsReal as calculated ...
const RooArgSet * nset() const
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
friend class RooRealSumFunc
RooNumIntConfig * specialIntegratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
RooAbsMoment * sigma(RooRealVar &obs)
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets,...
RooArgList is a container object that can hold multiple RooAbsArg objects.
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Class RooBinning is an implements RooAbsBinning in terms of an array of boundary values,...
double * array() const override
Return array of boundary values.
Int_t numBoundaries() const override
Return the number boundaries.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Getter function without integration set.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Setter function without integration set.
bool setLabel(const char *label, bool printError=true) override
Set value by specifying the name of the desired state.
RooChangeTracker is a meta object that tracks value changes in a given set of RooAbsArgs by registeri...
bool hasChanged(bool clearState)
Returns true if state has changed since last call with clearState=true.
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
RooConstVar represent a constant real-valued object.
RooCustomizer is a factory class to produce clones of a prototype composite PDF object with the same ...
void replaceArg(const RooAbsArg &orig, const RooAbsArg &subst)
Replace any occurence of arg 'orig' with arg 'subst'.
RooAbsArg * build(const char *masterCatState, bool verbose=false)
Build a clone of the prototype executing all registered 'replace' rules and 'split' rules for the mas...
RooLinearVar is the most general form of a derived real-valued object that can be used by RooRealInte...
void calculateFractions(const RooMomentMorphFuncND &self, bool verbose=true) const
RooChangeTracker * _tracker
RooArgList containedArgs(Action) override
std::vector< int > _nnuis
void addBinning(const RooAbsBinning &binning)
std::vector< RooAbsBinning * > _grid
std::vector< std::vector< double > > _nref
void addPdf(const RooAbsReal &func, int bin_x)
RooAbsReal * sumFunc(const RooArgSet *nset)
RooObjCacheManager _cacheMgr
! Transient cache manager
void findShape(const std::vector< double > &x) const
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
void initializeParameters(const RooArgList &parList)
static bool next_combination(const Iterator first, Iterator k, const Iterator last)
CacheElem * getCache(const RooArgSet *nset) const
static void cartesian_product(std::vector< std::vector< T > > &out, std::vector< std::vector< T > > &in)
void initializeObservables(const RooArgList &obsList)
bool setBinIntegrator(RooArgSet &allVars)
std::vector< std::vector< double > > _squareVec
virtual double getVal(const RooArgSet *set=0) const
int sij(const int &i, const int &j) const
~RooMomentMorphFuncND() override
TIterator * _obsItr
do not persist
std::vector< int > _squareIdx
const RooArgSet & getConfigSection(const char *name) const
Retrieve configuration information specific to integrator with given name.
RooRealVar represents a variable that can be changed from the outside.
virtual TObject * Next()=0
TMatrixT< Element > & Invert(Double_t *det=0)
Invert the matrix and calculate its determinant.
const char * GetName() const override
Returns name of object.
Element * GetMatrixArray()
void(off) SmallVectorTemplateBase< T
std::shared_ptr< std::function< double(double)> > Linear
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
std::vector< T >::const_iterator begin