Tutorial

After importing Gepard package, if one of the theories available in gepard.fits are not enough, because, for example, we want to make a new GPD fit, we can construct a new “theory” (more precissely, a subclass of the Theory class). We do this by combining code for GPDs, CFFs and observables like in the next example:

1. for GPDs we take SO(3) partial-waves conformal space model PWNormGPD,

2. for CFFs we take MellinBarnesCFF which combines conformal space GPDs with appropriate hard-scattering coefficients,

3. for implementation of DVCS observables we take BMK

4. we also use DIS to be able to check that forward limit of GPDs describes DIS data

This is combined like this:

>>> import gepard as g
>>> class MyTheory(g.PWNormGPD, g.MellinBarnesCFF, g.DIS, g.BMK):
...     pass

(Source code)

We then construct an instance of this new class, where we decide that we want to work at NLO (p = 1):

>>> th = MyTheory(p=1)

(Source code)

Note

It is somewhat of a convention, which we follow in this documentation, that names of Theory instances always start with th, while the names of DataPoint instances always start with pt.

We would now like to confront this theory with experimental data on DIS and DVCS. For demonstration purposes, we take just two H1 datasets:

>>> DISpoints = g.dset[206]
>>> DVCSpoints = g.dset[39]

(Source code)

To see in more detail what are these datasets, we can use utility function describe_data:

>>> g.describe_data(DISpoints+DVCSpoints)
npt x obs     collab  FTn    id  ref.
----------------------------------------------
10 x DISF2   H1      N/A    206 Nucl.Phys.B470(96)3
 8 x XGAMMA  H1      N/A    39  hep-ex/0505061
----------------------------------------------
TOTAL = 18

(Source code)

Here we see that dataset id=206 contains 10 measurements of DIS \(F_2\), while id=39 contains 8 measurements of pure DVCS cross-section, differential in \(t\). (FTn column is not relevant here.). We also see the reference to original literature where measurements were published.

To see if our theory describes this data we can either calculate total \(\chi^2\):

>>> th.chisq(DVCSpoints)
85.21

(Source code)

or we can plot theory line against data points like this:

>>> import gepard.plots
>>> gepard.plots.jbod(points=DVCSpoints, lines=th).show()

(Source code, png, hires.png, pdf)

This is obviously bad, so let us fit the parameters of the theory to this data. For this, we construct the MinuitFitter object, release some of the model parameters (overal normalization ns, residual \(t\)-dependence parameter ms2, and normalization of the second partial wave secs, all for sea quarks):

>>> f = g.MinuitFitter(DISpoints+DVCSpoints, th)
>>> f.release_parameters('ns', 'ms2', 'secs')
>>> f.fit()

(Source code)

After fitting is done, we print the resulting values and uncertainties of fitting parameters:

>>> th.print_parameters()
ns    =    0.17 +- 0.01
ms2   =    0.93 +- 0.10
secs  =    0.18 +- 0.03

(Source code)

Theory now describes the data fine, as one can see from \(\chi^2\) value:

>>> th.chisq(DISpoints+DVCSpoints)
6.33

(Source code)

and, visually, from the plot:

>>> gepard.plots.jbod(points=DVCSpoints, lines=th).show()

(Source code, png, hires.png, pdf)

Finally, one could calculate and then plot some particular CFF, like this:

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> xis = np.linspace(0.001, 0.1)
>>> ims = []
>>> res = []
>>> for xi in xis:
...     pt = g.DataPoint(xi=xi, t=-0.2, Q2=4)
...     ims.append(xi*th.ImH(pt))
...     res.append(xi*th.ReH(pt))
>>> plt.plot(xis, ims, label='Im(H)')  
>>> plt.plot(xis, res, label='Re(H)')  
>>> plt.xlabel(r'$\xi$', fontsize=14)  
>>> plt.ylabel('ImH', fontsize=14)  
>>> plt.legend()  

(Source code, png, hires.png, pdf)