Logo ROOT  
Reference Guide
RooLagrangianMorphFunc.h
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * *
4 * authors: *
5 * Lydia Brenner (lbrenner@cern.ch), Carsten Burgard (cburgard@cern.ch) *
6 * Katharina Ecker (kecker@cern.ch), Adam Kaluza (akaluza@cern.ch) *
7 * Copyright (c) 2000-2007, Regents of the University of California *
8 * and Stanford University. All rights reserved. *
9 * *
10 * Redistribution and use in source and binary forms, *
11 * with or without modification, are permitted according to the terms *
12 * listed in LICENSE (http://roofit.sourceforge.net/license.txt)
13 *****************************************************************************/
14
15////////////////////////////////////////////////////////////////////////////////////////////////
16//
17// RooLagrangianMorphFunc
18//
19// The RooLagrangianMorphFunc is a type of RooAbsReal that allows to morph
20// different input EFT samples to some arbitrary output EFT
21// sample, as long as the desired set of output parameters lie
22// within the realm spanned by the input samples. More
23// specifically, it expects as an input a TFile (or TDirectory)
24// with the following layout:
25//
26// TDirectory
27// |-sample1
28// | |-param_card // TH1 EFT parameter values of sample1
29// | | histogram1 // TH1 of some physics distribution
30// | |-subfolder1 // a subfolder (optional)
31// | | |-histogram2 // TH1 of some physics distribution
32// | | |-....
33// |-sample2
34// | |-param_card // TH1 of EFT parameter values of sample2
35// | | histogram1 // TH1 of some physics distribution
36// | |-subfolder1 // same folder structure as before
37// | | |-histogram2 // TH1 of some physics distribution
38// | | |-....
39// |-sampleN
40// The RooLagrangianMorphFunc operates on this structure, extracts data
41// and meta-data and produces a morphing result as a RooRealSumFunc
42// consisting of the input histograms with appropriate prefactors.
43//
44// The histograms to be morphed can be accessed via their paths in
45// the respective sample, e.g. using
46// "histogram"
47// or "subfolder1/histogram1"
48// or "some/deep/path/to/some/subfolder/histname"
49//
50////////////////////////////////////////////////////////////////////////////////////////////////
51
52#ifndef ROO_LAGRANGIAN_MORPH
53#define ROO_LAGRANGIAN_MORPH
54
55#include "RooFit/Floats.h"
56#include "RooAbsArg.h"
57#include "RooAbsReal.h"
58#include "RooArgList.h"
59#include "RooRatio.h"
60#include "RooRealSumFunc.h"
61#include "RooRealSumPdf.h"
62#include "RooSetProxy.h"
63#include "RooWrapperPdf.h"
64#include "TMatrixD.h"
65
66class RooWorkspace;
68class RooProduct;
69class RooRealVar;
70class TPair;
71class TFolder;
73
74#include <fstream>
75#include <iostream>
76#include <string>
77#include <vector>
78
80
81public:
82 typedef std::map<const std::string, double> ParamSet;
83 typedef std::map<const std::string, int> FlagSet;
84 typedef std::map<const std::string, ParamSet> ParamMap;
85 typedef std::map<const std::string, FlagSet> FlagMap;
86
87 struct Config {
88
89 std::string observableName;
90 std::string fileName;
93 std::vector<std::string> folderNames;
98 std::vector<RooArgList *> vertices;
99 std::vector<RooArgList *> nonInterfering;
101 };
102
104 RooLagrangianMorphFunc(const char *name, const char *title, const char *filename, const char *observableName,
105 const RooArgSet &couplings, const RooArgSet &inputs);
106 RooLagrangianMorphFunc(const char *name, const char *title, const Config &config);
107 RooLagrangianMorphFunc(const RooLagrangianMorphFunc &other, const char *newName);
108
109 ~RooLagrangianMorphFunc() override;
110
111 std::list<double> *
112 binBoundaries(RooAbsRealLValue & /*obs*/, double /*xlo*/, double /*xhi*/) const override;
113 std::list<double> *
114 plotSamplingHint(RooAbsRealLValue & /*obs*/, double /*xlo*/, double /*xhi*/) const override;
115 bool isBinnedDistribution(const RooArgSet &obs) const override;
116 double evaluate() const override;
117 TObject *clone(const char *newname) const override;
118 double getValV(const RooArgSet *set = 0) const override;
119
120 bool checkObservables(const RooArgSet *nset) const override;
121 bool forceAnalyticalInt(const RooAbsArg &arg) const override;
122 Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet,
123 const char *rangeName = 0) const override;
124 double
125 analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName = 0) const override;
126 void printMetaArgs(std::ostream &os) const override;
127 RooAbsArg::CacheMode canNodeBeCached() const override;
128 void setCacheAndTrackHints(RooArgSet &) override;
129
131
132 void setParameters(const char *foldername);
133 void setParameters(TH1 *paramhist);
134 void setParameter(const char *name, double value);
135 void setFlag(const char *name, double value);
136 void setParameters(const ParamSet &params);
137 void setParameters(const RooArgList *list);
138 double getParameterValue(const char *name) const;
139 RooRealVar *getParameter(const char *name) const;
140 RooRealVar *getFlag(const char *name) const;
141 bool hasParameter(const char *paramname) const;
142 bool isParameterUsed(const char *paramname) const;
143 bool isParameterConstant(const char *paramname) const;
144 void setParameterConstant(const char *paramname, bool constant) const;
145 void setParameter(const char *name, double value, double min, double max);
146 void setParameter(const char *name, double value, double min, double max, double error);
147 void randomizeParameters(double z);
148 const RooArgSet *getParameterSet() const;
149 ParamSet getMorphParameters(const char *foldername) const;
151
153
154 int nParameters() const;
155 int nPolynomials() const;
156
157 bool isCouplingUsed(const char *couplname);
158 const RooArgList *getCouplingSet() const;
159 ParamSet getCouplings() const;
160
161 TMatrixD getMatrix() const;
163 double getCondition() const;
164
165 RooRealVar *getObservable() const;
166 RooRealVar *getBinWidth() const;
167
168 void printEvaluation() const;
169 void printCouplings() const;
170 void printFlags() const;
171 void printPhysics() const;
172
173 RooProduct *getSumElement(const char *name) const;
174
175 std::vector<std::string> getSamples() const;
176
177 double expectedUncertainty() const;
178 TH1 *createTH1(const std::string &name);
179 TH1 *createTH1(const std::string &name, bool correlateErrors);
180
181protected:
182 class CacheElem;
183 void init();
184 void setup(bool ownParams = true);
185 bool _ownParameters = false;
186
187 mutable RooObjCacheManager _cacheMgr; //! The cache manager
188
189 void addFolders(const RooArgList &folders);
190
191 bool hasCache() const;
193 void updateSampleWeights();
194
195 RooRealVar *setupObservable(const char *obsname, TClass *mode, TObject *inputExample);
196
197public:
198 /// length of floating point digits precision supported by implementation.
199 static constexpr double implementedPrecision = RooFit::SuperFloatPrecision::digits10;
200
201 void writeMatrixToFile(const TMatrixD &matrix, const char *fname);
202 void writeMatrixToStream(const TMatrixD &matrix, std::ostream &stream);
203 TMatrixD readMatrixFromFile(const char *fname);
204 TMatrixD readMatrixFromStream(std::istream &stream);
205
206 int countSamples(std::vector<RooArgList *> &vertices);
207 int countSamples(int nprod, int ndec, int nboth);
208
209 TPair *makeCrosssectionContainer(double xs, double unc);
210 std::map<std::string, std::string>
211 createWeightStrings(const ParamMap &inputs, const std::vector<std::vector<std::string>> &vertices);
212 std::map<std::string, std::string>
213 createWeightStrings(const ParamMap &inputs, const std::vector<RooArgList *> &vertices, RooArgList &couplings);
214 std::map<std::string, std::string>
215 createWeightStrings(const ParamMap &inputs, const std::vector<RooArgList *> &vertices, RooArgList &couplings,
216 const FlagMap &flagValues, const RooArgList &flags,
217 const std::vector<RooArgList *> &nonInterfering);
218 RooArgSet createWeights(const ParamMap &inputs, const std::vector<RooArgList *> &vertices, RooArgList &couplings,
219 const FlagMap &inputFlags, const RooArgList &flags,
220 const std::vector<RooArgList *> &nonInterfering);
221 RooArgSet createWeights(const ParamMap &inputs, const std::vector<RooArgList *> &vertices, RooArgList &couplings);
222
223 bool updateCoefficients();
224 bool useCoefficients(const TMatrixD &inverse);
225 bool useCoefficients(const char *filename);
226 bool writeCoefficients(const char *filename);
227
228 int countContributingFormulas() const;
229 RooAbsReal *getSampleWeight(const char *name);
230 void printParameters(const char *samplename) const;
231 void printParameters() const;
232 void printSamples() const;
233 void printSampleWeights() const;
234 void printWeights() const;
235
236 void setScale(double val);
237 double getScale();
238
239 int nSamples() const { return this->_config.folderNames.size(); }
240
241 RooRealSumFunc *getFunc() const;
242 std::unique_ptr<RooWrapperPdf> createPdf() const;
243
245 double expectedEvents(const RooArgSet *nset) const;
246 double expectedEvents(const RooArgSet &nset) const;
247 double expectedEvents() const;
248 bool selfNormalized() const { return true; }
249
252
253 static std::unique_ptr<RooRatio> makeRatio(const char *name, const char *title, RooArgList &nr, RooArgList &dr);
254
255protected:
256 double _scale = 1;
257 std::map<std::string, int> _sampleMap;
264 std::vector<std::vector<RooListProxy *>> _diagrams;
265 mutable const RooArgSet *_curNormSet = nullptr; //!
266
267 // TODO: the _nonInterfering is not filled anywhere and also not considered
268 // correctly in the copy constructor. Can it be removed?
269 std::vector<RooListProxy *> _nonInterfering;
270
272};
273
274#endif
#define f(i)
Definition: RSha256.hxx:104
#define ClassDefOverride(name, id)
Definition: Rtypes.h:339
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 Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char filename
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char mode
char name[80]
Definition: TGX11.cxx:110
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:77
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
Class RooLagrangianMorphing is a implementation of the method of Effective Lagrangian Morphing,...
static constexpr double implementedPrecision
length of floating point digits precision supported by implementation.
bool isParameterConstant(const char *paramname) const
return true if the parameter with the given name is set constant, false otherwise
TPair * makeCrosssectionContainer(double xs, double unc)
create TPair containers of the type expected by the RooLagrangianMorph
bool isBinnedDistribution(const RooArgSet &obs) const override
check if this PDF is a binned distribution in the given observable
int nPolynomials() const
return the number of samples in this morphing function
void setParameter(const char *name, double value)
set one parameter to a specific value
ParamSet getMorphParameters() const
retrieve the parameter set
double evaluate() const override
call getVal on the internal function
RooProduct * getSumElement(const char *name) const
return the RooProduct that is the element of the RooRealSumPdfi corresponding to the given sample nam...
void insert(RooWorkspace *ws)
RooRealVar * getBinWidth() const
retrieve the histogram observable
void writeMatrixToFile(const TMatrixD &matrix, const char *fname)
write a matrix to a file
RooRealVar * getParameter(const char *name) const
retrieve the RooRealVar object incorporating the parameter with the given name
bool useCoefficients(const TMatrixD &inverse)
setup the morphing function with a predefined inverse matrix call this function before any other afte...
const RooArgSet * getParameterSet() const
get the set of parameters
std::map< const std::string, int > FlagSet
TMatrixD readMatrixFromStream(std::istream &stream)
read a matrix from a stream
std::vector< std::string > getSamples() const
return the vector of sample names, used to build the morph func
int countSamples(std::vector< RooArgList * > &vertices)
calculate the number of samples needed to morph a certain physics process usage: countSamples ( { Roo...
void setCacheAndTrackHints(RooArgSet &) override
Retrieve the matrix of coefficients.
ParamSet getCouplings() const
retrieve a set of couplings (-?-)
void printSampleWeights() const
print the current sample weights
std::map< const std::string, double > ParamSet
RooArgSet createWeights(const ParamMap &inputs, const std::vector< RooArgList * > &vertices, RooArgList &couplings, const FlagMap &inputFlags, const RooArgList &flags, const std::vector< RooArgList * > &nonInterfering)
create only the weight formulas. static function for external usage.
RooLagrangianMorphFunc * getLinear() const
void writeMatrixToStream(const TMatrixD &matrix, std::ostream &stream)
write a matrix to a stream
void addFolders(const RooArgList &folders)
The cache manager.
std::map< const std::string, ParamSet > ParamMap
bool updateCoefficients()
retrive the new physics objects and update the weights in the morphing function
RooRealVar * getObservable() const
retrieve the histogram observable
int countContributingFormulas() const
count the number of formulas that correspond to the current parameter set
std::list< double > * plotSamplingHint(RooAbsRealLValue &, double, double) const override
retrieve the sample Hint
RooRealVar * getFlag(const char *name) const
retrieve the RooRealVar object incorporating the flag with the given name
void randomizeParameters(double z)
randomize the parameters a bit useful to test and debug fitting
bool isCouplingUsed(const char *couplname)
check if there is any morphing power provided for the given coupling morphing power is provided as so...
void readParameters(TDirectory *f)
read the parameters from the input file
double getScale()
get energy scale of the EFT expansion
std::vector< RooListProxy * > _nonInterfering
double getCondition() const
Retrieve the condition of the coefficient matrix.
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const override
Retrieve the matrix of coefficients.
TMatrixD getMatrix() const
Retrieve the matrix of coefficients.
void printWeights() const
print the current sample weights
void printCouplings() const
print a set of couplings
TMatrixD readMatrixFromFile(const char *fname)
read a matrix from a text file
~RooLagrangianMorphFunc() override
default destructor
void printParameters() const
print the parameters and their current values
void printPhysics() const
print the current physics values
static std::unique_ptr< RooRatio > makeRatio(const char *name, const char *title, RooArgList &nr, RooArgList &dr)
Return the RooRatio form of products and denominators of morphing functions.
void setFlag(const char *name, double value)
set one flag to a specific value
TH1 * createTH1(const std::string &name)
retrieve a histogram output of the current morphing settings
double expectedUncertainty() const
return the expected uncertainty for the current parameter set
int nParameters() const
return the number of parameters in this morphing function
RooLagrangianMorphFunc::CacheElem * getCache(const RooArgSet *nset) const
retrieve the cache object
bool hasParameter(const char *paramname) const
check if a parameter of the given name is contained in the list of known parameters
bool checkObservables(const RooArgSet *nset) const override
check if observable exists in the RooArgSet (-?-)
bool hasCache() const
return true if a cache object is present, false otherwise
void printFlags() const
print the flags and their current values
void setScale(double val)
set energy scale of the EFT expansion
TMatrixD getInvertedMatrix() const
Retrieve the matrix of coefficients after inversion.
RooAbsArg::CacheMode canNodeBeCached() const override
Retrieve the matrix of coefficients.
void updateSampleWeights()
update sample weight (-?-)
void setParameters(const char *foldername)
set the morphing parameters to those supplied in the sample with the given name
bool isParameterUsed(const char *paramname) const
check if there is any morphing power provided for the given parameter morphing power is provided as s...
RooAbsPdf::ExtendMode extendMode() const
return extended mored capabilities
std::map< const std::string, FlagSet > FlagMap
bool forceAnalyticalInt(const RooAbsArg &arg) const override
Force analytical integration for the given observable.
void setParameterConstant(const char *paramname, bool constant) const
call setConstant with the boolean argument provided on the parameter with the given name
void printSamples() const
print all the known samples to the console
double getParameterValue(const char *name) const
set one parameter to a specific value
double getValV(const RooArgSet *set=0) const override
call getVal on the internal function
void setup(bool ownParams=true)
setup this instance with the given set of operators and vertices if own=true, the class will own the ...
void printMetaArgs(std::ostream &os) const override
Retrieve the matrix of coefficients.
std::unique_ptr< RooWrapperPdf > createPdf() const
(currently similar to cloning the Pdf
RooAbsReal * getSampleWeight(const char *name)
retrieve the weight (prefactor) of a sample with the given name
std::map< std::string, std::string > createWeightStrings(const ParamMap &inputs, const std::vector< std::vector< std::string > > &vertices)
create only the weight formulas. static function for external usage.
std::vector< std::vector< RooListProxy * > > _diagrams
double expectedEvents() const
return the number of expected events for the current parameter set
std::map< std::string, int > _sampleMap
TObject * clone(const char *newname) const override
cloning method
RooRealVar * setupObservable(const char *obsname, TClass *mode, TObject *inputExample)
setup observable, recycle existing observable if defined
const RooArgList * getCouplingSet() const
get the set of couplings
RooRealSumFunc * getFunc() const
get the func
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
retrieve the list of bin boundaries
void printEvaluation() const
print the contributing samples and their respective weights
bool writeCoefficients(const char *filename)
write the inverse matrix to a file
void collectInputs(TDirectory *f)
retrieve the physics inputs
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const override
Retrieve the mat.
Class RooObjCacheManager is an implementation of class RooCacheManager<RooAbsCacheElement> and specia...
A histogram function that assigns scale parameters to every bin.
A RooProduct represents the product of a given set of RooAbsReal objects.
Definition: RooProduct.h:29
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:43
TClass instances represent classes, structs and namespaces in the ROOT type system.
Definition: TClass.h:80
Describe directory structure in memory.
Definition: TDirectory.h:45
A TFolder object is a collection of objects and folders.
Definition: TFolder.h:30
TH1 is the base class of all histogram classes in ROOT.
Definition: TH1.h:58
Mother of all ROOT objects.
Definition: TObject.h:37
Class used by TMap to store (key,value) pairs.
Definition: TMap.h:102
std::vector< RooArgList * > nonInterfering
std::vector< std::string > folderNames
std::vector< RooArgList * > vertices
void ws()
Definition: ws.C:66