xsec interpolation

[tomeichlersmith]

The actual interpolation method is rather simple - we choose the three nearest points and then use a quadratic interpolation scheme. This makes the interpolation extremely fast.

Following Michael R's advice, I have the initial table of energy/cross section values be determined upon the first "table-miss" so there isn't much of a speed up in low-event-number runs; however, constructing this table immediately allows for a majority of following cross section calls to "hit" the interpolation instead of the full calculation, saving a lot of time.

My validation is copied below, showing that the interpolation is accurate. I still need to integrate this into the simulation so that I can check its speed.

To Do

• [x] Integrate into simulation
• [x] Check speed of interpolation vs cache
• [x] Write dump method so people can inspect the sample points being used in the interpolation

---

import pandas

%matplotlib inline
import matplotlib.pyplot as plt
import mplhep
plt.style.use(mplhep.style.ROOT)

Validation of Interpolation

I want to compare the full cross section calculation to the interpolation. For this to function, I do the full calculation at a set of points and then run the interpolation through those points at a higher granularity so that we see the interpolation do its "magic".

After building this branch, I ran the following bash script to calculate the necessary CSV data files.

for method in fullww iww hiww; do
  for mass in 1 0.1 0.01 0.001; do
    ./g4db-xsec-calc --apmass ${mass} --method ${method} --energy 0 4 0.1 --output no-interpo-${mass}GeV-${method}-xsec.csv
    ./g4db-xsec-calc --apmass ${mass} --method ${method} --energy 0 4 0.01 --interpolation --output interpo-${mass}GeV-${method}-xsec.csv
  done
done

for mass in [1,0.1,0.01,0.001] :
    (raw, ratio) = plt.gcf().subplots(
        ncols = 1, nrows = 2, 
        sharex = 'col', gridspec_kw=dict(height_ratios = [3,1])
    )
    plt.subplots_adjust(hspace=0)

    for method in ['full','i','hi'] :
        no_interpo = pandas.read_csv(f'no-interpo-{method}ww-{mass}GeV-xsec.csv')
        interpo    = pandas.read_csv(f'interpo-{method}ww-{mass}GeV-xsec.csv')
        raw.plot(interpo['Energy [MeV]'], interpo['Xsec [pb]'], label=f'{method} interpo')
        raw.scatter(no_interpo['Energy [MeV]'], no_interpo['Xsec [pb]'], label=f'{method} nointerpo')
        # since the interpolated data set and the non-interpolated data sets have different data points
        #   we use the non-interpolated energies to try to find a selection for the interpolated ones
        # also since we are dividing by no_interpo, we select only those energies with a positive xsec
        selection = interpo['Energy [MeV]'].isin(no_interpo[no_interpo['Xsec [pb]']>0]['Energy [MeV]'])
        ratio.plot(
            interpo[selection]['Energy [MeV]'], 
            interpo[selection]['Xsec [pb]'].to_numpy()/no_interpo[no_interpo['Xsec [pb]']>0]['Xsec [pb]'].to_numpy()
        )
    raw.legend(
        title=f'$m_A = {mass}$ GeV'
    )
    raw.set_ylabel('Xsec [pb]')
    ratio.set_ylabel('Interp Acc')
    ratio.set_xlabel('Energy [MeV]')
    plt.show()