Interactive fitting with Iminuit¶

This example uses Moose in combination with the iminuit (DOI) fitting library, maintained by CERN's ROOT team.

Iminuit provides its own interface for interactive plotting, relying on matplotlib.

It works however similarly to those in the other interactive examples, though does not update the plot as frequently.

The fit procedure can be performed on either of two included sample files in the folder ./N2SPS:

1. A spectrum of the second positive system (0,0) band, as recorded at the start of a nanosecond repetitively pulsed discharge: 20220302_s01_sig.cam

2. A spectrum simulated by SpecAir as part of the unit test from RADIS

   • See also: the corresponding GitHub page
   • Or the notebook by the RADIS developers that demonstrates using MassiveOES

In [1]:

from pathlib import Path
import pandas as pd
import matplotlib.pyplot as plt
from iminuit import Minuit
from iminuit.cost import LeastSquares

import Moose
from Moose.iminuit import as_minuit_objective

db = Moose.query_DB('N2CB')
fs =list(Path("N2SPS").glob('*'))

from pathlib import Path import pandas as pd import matplotlib.pyplot as plt from iminuit import Minuit from iminuit.cost import LeastSquares import Moose from Moose.iminuit import as_minuit_objective db = Moose.query_DB('N2CB') fs =list(Path("N2SPS").glob('*'))

Preliminary setup¶

To work with iminuit we need to change the standard Moose.model_for_fit function signature, such that it conforms to what iminuit expects.

The most inconvenient aspect is that iminuit expects all keyword arguments to be floats, which is not the case.

To achieve this, you can use the Moose.iminuit.minuit_objective_decorator decorator/wrapper, to inject the non-float kwargs into the function before iminuit inspects it.

Additionally, for clarity, the plot_func function is defined below, to override the default way that iminuit uses to plot data during the interactive process.

In [2]:

# Define the objective function
objective = as_minuit_objective(Moose.model_for_fit,sim_db = db, normalize=True)

def plot_func(params,x=None,y=None,show_params=False, ):
    """A custom plotting function to override the default way iminuit chooses to plot data/fit.

    The `objective_fn` argument must be the objective function to be fitted, so that it can be called to plot the fit result.
    
    Make sure to call the `interactive` method of the Minuit object as: `m.interactive(plot=plot_func,x=...,y=..., objective_fn=...)`
    """
    plt.gca().plot(x,y, label="Data")
    plt.gca().plot(x,objective(x,*params), label='Fit')
    plt.gca().legend()
    if show_params:
        for i,(name, value) in enumerate(params.to_dict().items()):
            plt.gca().annotate(f"{name}:{value:.3g}",(0.05,0.92-i*0.08), xycoords='axes fraction')

# Define the objective function objective = as_minuit_objective(Moose.model_for_fit,sim_db = db, normalize=True) def plot_func(params,x=None,y=None,show_params=False, ): """A custom plotting function to override the default way iminuit chooses to plot data/fit. The `objective_fn` argument must be the objective function to be fitted, so that it can be called to plot the fit result. Make sure to call the `interactive` method of the Minuit object as: `m.interactive(plot=plot_func,x=...,y=..., objective_fn=...)` """ plt.gca().plot(x,y, label="Data") plt.gca().plot(x,objective(x,*params), label='Fit') plt.gca().legend() if show_params: for i,(name, value) in enumerate(params.to_dict().items()): plt.gca().annotate(f"{name}:{value:.3g}",(0.05,0.92-i*0.08), xycoords='axes fraction')

In [3]:

# Load the data
sample_1 = pd.read_csv(fs[0], sep='\t', header=1, names=['Wavelength', 'I'])
sample_1['Norm'] = (sample_1['I']-sample_1['I'].min())/(sample_1['I'].max()-sample_1['I'].min())

# Setup iminuit with model and apply parameter bounds/fixing
ls_1 = LeastSquares(sample_1.iloc[:,0],sample_1.iloc[:,2],sample_1.iloc[:,2].nsmallest(100).std(), objective)
m_1 = Minuit(ls_1, **{k:v['value'] for k,v in Moose.default_params.items()})
for k,v in Moose.default_params.items():
    if v.get('vary',True):
        m_1.limits[k] = (v['min'], v['max'])
    else:
        m_1.fixto(k,v['value'])

# Visualize an interactive interface for fitting
m_1.interactive(plot=plot_func, x=sample_1.iloc[:,0],y = sample_1.iloc[:,2])

# Load the data sample_1 = pd.read_csv(fs[0], sep='\t', header=1, names=['Wavelength', 'I']) sample_1['Norm'] = (sample_1['I']-sample_1['I'].min())/(sample_1['I'].max()-sample_1['I'].min()) # Setup iminuit with model and apply parameter bounds/fixing ls_1 = LeastSquares(sample_1.iloc[:,0],sample_1.iloc[:,2],sample_1.iloc[:,2].nsmallest(100).std(), objective) m_1 = Minuit(ls_1, **{k:v['value'] for k,v in Moose.default_params.items()}) for k,v in Moose.default_params.items(): if v.get('vary',True): m_1.limits[k] = (v['min'], v['max']) else: m_1.fixto(k,v['value']) # Visualize an interactive interface for fitting m_1.interactive(plot=plot_func, x=sample_1.iloc[:,0],y = sample_1.iloc[:,2])

Out[3]:

In [4]:

# Load the data
sample_2 = pd.read_csv(fs[1], sep='\t', header=1, names=['Wavelength', 'I'])
sample_2['Norm'] = (sample_2['I']-sample_2['I'].min())/(sample_2['I'].max()-sample_2['I'].min())

# Setup iminuit with model and apply parameter bounds/fixing
ls_2 = LeastSquares(sample_2.iloc[:,0],sample_2.iloc[:,2],sample_2.iloc[:,2].nsmallest(100).std(), objective)
m_2 = Minuit(ls_2, **{k:v['value'] for k,v in Moose.default_params.items()})

for k,v in Moose.default_params.items():
    if v.get('vary',True):
        m_2.limits[k] = (v['min'], v['max'])
    else:
        m_2.fixto(k,v['value'])

# Visualize an interactive interface for fitting, now showing the parameter values
m_2.interactive(plot=plot_func,x= sample_2.iloc[:,0], y=sample_2.iloc[:,2], show_params=True)

# Load the data sample_2 = pd.read_csv(fs[1], sep='\t', header=1, names=['Wavelength', 'I']) sample_2['Norm'] = (sample_2['I']-sample_2['I'].min())/(sample_2['I'].max()-sample_2['I'].min()) # Setup iminuit with model and apply parameter bounds/fixing ls_2 = LeastSquares(sample_2.iloc[:,0],sample_2.iloc[:,2],sample_2.iloc[:,2].nsmallest(100).std(), objective) m_2 = Minuit(ls_2, **{k:v['value'] for k,v in Moose.default_params.items()}) for k,v in Moose.default_params.items(): if v.get('vary',True): m_2.limits[k] = (v['min'], v['max']) else: m_2.fixto(k,v['value']) # Visualize an interactive interface for fitting, now showing the parameter values m_2.interactive(plot=plot_func,x= sample_2.iloc[:,0], y=sample_2.iloc[:,2], show_params=True)

Out[4]:

In [ ]: