RooFFTConvPdf: Resolution depends on convolution observable

Roofit and RooStats
1

Hi,

I have a question regarding the behaviour of RooFFTConvPdf which probably stems from an inherent misunderstanding of how the convolution process works, but here goes.

I want to convolve a “physics” pdf with another “resolution” pdf which depends on the value of the physics pdf observable. More concretely, if the observable is x, I want to convolve with a gaussian with a width which varies as a function of x. i.e.

f(x) (x) g(x’;mean = 0, sigma(x))

Here x and x’ refer to the same observable i.e. the same RooRealVar, but I don’t think they should be interpreted as identical for sampling of the resolution pdf. Accordingly, my question relates to how RooFFTConvPdf samples the input pdfs. Will it resample the resolution pdf for every value of x and correctly substitute that value of x into sigma (i.e. effectively convolving with different symmetric gaussian at each x) or will it instead sample the resolution pdf once and calculate sigma(x’) (i.e. convolving with the same skewed gaussian at all x)? Here’s come code which may do what I want (it compiles and runs fine at least!).

[code]// Set Up Common Variables
RooRealVar energy(“energy”, “Photon Energy”, 1e-3, 200e-3, “TeV”);

// Set up Power law physics model
RooRealVar plIndex(“plIndex”, “Power Law Spectral Index”, -2, -100000, -0.01);
RooGenericPdf spec(“spec”, “pow(@0,@1)”, RooArgList(energy, plIndex));

// Set up gaussian resolution model
RooProduct resSigNC(“resSigNC”, “Resolution Model Sigma Without Offset”, RooArgSet(energy, RooConst(0.1)));
RooAddition resSigma(“resSigma”, “Resolution Model Sigma”, RooArgSet(resSigNC, RooConst(1e-4)));

RooGaussian resMod(“resMod”, “Resolution Model”, energy, RooConst(0.0), resSigma);

// Set energy binning to adequate resolution
energy.setBins(10000, “cache”);

// Set up smeared power law model
RooFFTConvPdf smearedSpec(“smearedSpec”, “Smeared Spectrum”, energy, spec, resMod);
smearedSpec.setBufferFraction(1.0);
smearedSpec.setBufferStrategy(RooFFTConvPdf::Flat);

RooPlot * eFrame = energy.frame();
spec.plotOn(eFrame);
smearedSpec.plotOn(eFrame, LineColor(kRed));[/code]

It’s not completely clear to me whether RooFFTConvPdf is doing what I want or not. Perhaps what I want to do is in fact a 2D convolution?

Any help would be much appreciated,

Regards,

Hugh

2

Hi Hugh,

The FFT convolution pdf does not handle the case you describe correctly. It it samples
the pdf once (for a given set of parameter values) and then performes the convolution.
I’m not sure if I can make it work any differently because because of the way the
FFT procedure works.

Wouter

3

Hi Wouter,

Thanks for the reply. I thought this might be the case. Never mind.

Cheers,

Hugh

4

I have exactly the same problem. Are you trying to convolve energy with a resolution depending on the energy?

I solved with:

  1. create a bidimensional model
  2. project

see attachment
roofit_prova.C (2.1 KB)

5

Hi,

Thanks for the advice. In fact I already tried this approach - it’s similar to RooFit’s built in “brute force” numerical integration method. However, my problem was that I needed to perform this operation many thousands of times for a toy MC study and numerically integrating a power law just took too long. In the end, I decided to discretize the problem, and use a precomputed migration matrix.

Thanks again,

Hugh

6

Hi,

I have a question for you wiso.

I have slightly modified your script in order to fit the projected model “result” to the data “reco_data”.
result->fitTo(*reco_data);

It seems that the observable mtrue is considered as a variable in the fit whereas it should not.

[size=85]COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=8959.65 FROM HESSE STATUS=OK 12 CALLS 110 TOTAL
EDM=0.33296 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 mtrue 3.00000e+04 5.15036e+02 1.16430e-04 2.44171e+03
2 true_width 5.15049e+02 3.34128e+01 8.80727e-05 -1.36746e+00[/size]

How could I indicate to RooFit that mtrue is actually not a parameter of the model to fit ?
Shuld I set mtrue->SetConstant() ? (altough I think it is a bit brute).

Cheers,

Mathieu

7

HI Matthieu,

It looks to me that true is a conditional observables. You would need an input data set in true in order to fit your data and then flag it as ConditionalObservable(true) in the fitTo command.
See for example the conditional fit examples in roofit/tutorial
root.cern.ch/root/htmldoc/tutor … ors.C.html

Best Regards

Lorenzo

8

Hi Lorenzo,
thanks for the reply.
The point is that I only have a dataset in reco and wants to fit this data with a folded model.
What I do not understand is why Roofit interprets true as a conditional as the model I want to fit is a projection over true (= integration over true).

Cheers,
Mathieu

9

Hi,

I am not sure I have fully understood your model. If “true” is a variable which is integrated over, then you should not have any more dependency on that variable.
If the model somehow depends still on it, you need either to have a data set or a distribution (a pdf ) on that variable. In the first case is a conditional variable, in the second case, it can be either a parameter of the model or it could be a “global observable”. If it is a global observable, you need to set it constant when fitting to its nominal value.

Lorenzo