Skip to content

Commit

Permalink
Use iminuit>=2
Browse files Browse the repository at this point in the history
  • Loading branch information
nbiederbeck committed May 11, 2021
1 parent b4d1e29 commit 524e621
Show file tree
Hide file tree
Showing 4 changed files with 31 additions and 42 deletions.
26 changes: 11 additions & 15 deletions ctapipe/image/muon/fitting.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,7 @@

from ctapipe.utils.quantities import all_to_value

__all__ = [
"kundu_chaudhuri_circle_fit",
"taubin_circle_fit",
]
__all__ = ["kundu_chaudhuri_circle_fit", "taubin_circle_fit"]


def kundu_chaudhuri_circle_fit(x, y, weights):
Expand Down Expand Up @@ -79,18 +76,17 @@ def taubin_circle_fit(x, y, mask):

# minimization method
fit = Minuit(
make_taubin_loss_function(x_masked, y_masked),
xc=xc,
yc=yc,
r=taubin_r_initial,
error_xc=taubin_error,
error_yc=taubin_error,
error_r=taubin_error,
limit_xc=(-2 * R, 2 * R),
limit_yc=(-2 * R, 2 * R),
limit_r=(0, R),
pedantic=False,
make_taubin_loss_function(x_masked, y_masked), xc=xc, yc=yc, r=taubin_r_initial
)

fit.errors["xc"] = taubin_error
fit.errors["yc"] = taubin_error
fit.errors["r"] = taubin_error

fit.limits["xc"] = (-2 * R, 2 * R)
fit.limits["yc"] = (-2 * R, 2 * R)
fit.limits["r"] = (0, R)

fit.migrad()

radius = Quantity(fit.values["r"], original_unit)
Expand Down
43 changes: 18 additions & 25 deletions ctapipe/image/muon/intensity_fitter.py
Original file line number Diff line number Diff line change
Expand Up @@ -517,37 +517,30 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non
pedestal=pedestal,
hole_radius=self.hole_radius_m.tel[tel_id] * u.m,
)
negative_log_likelihood.errordef = Minuit.LIKELIHOOD

initial_guess = create_initial_guess(center_x, center_y, radius, telescope)

step_sizes = {}
step_sizes["error_impact_parameter"] = 0.5
step_sizes["error_phi"] = np.deg2rad(0.5)
step_sizes["error_ring_width"] = 0.001 * radius.to_value(u.rad)
step_sizes["error_optical_efficiency_muon"] = 0.05

constraints = {}
constraints["limit_impact_parameter"] = (0, None)
constraints["limit_phi"] = (-np.pi, np.pi)
constraints["fix_radius"] = True
constraints["fix_center_x"] = True
constraints["fix_center_y"] = True
constraints["limit_ring_width"] = (0.0, None)
constraints["limit_optical_efficiency_muon"] = (0.0, None)

# Create Minuit object with first guesses at parameters
# strip away the units as Minuit doesnt like them

minuit = Minuit(
negative_log_likelihood,
# forced_parameters=parameter_names,
**initial_guess,
**step_sizes,
**constraints,
errordef=0.5,
print_level=0,
pedantic=True,
)
minuit = Minuit(negative_log_likelihood, **initial_guess)

minuit.print_level = 0

minuit.errors["impact_parameter"] = 0.5
minuit.errors["phi"] = np.deg2rad(0.5)
minuit.errors["ring_width"] = 0.001 * radius.to_value(u.rad)
minuit.errors["optical_efficiency_muon"] = 0.05

minuit.limits["impact_parameter"] = (0, None)
minuit.limits["phi"] = (-np.pi, np.pi)
minuit.limits["ring_width"] = (0.0, None)
minuit.limits["optical_efficiency_muon"] = (0.0, None)

minuit.fixed["radius"] = True
minuit.fixed["center_x"] = True
minuit.fixed["center_y"] = True

# Perform minimisation
minuit.migrad()
Expand Down
2 changes: 1 addition & 1 deletion environment.yml
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ dependencies:
- cython
- graphviz
- h5py
- iminuit>=1.3,<2.0.0a0
- iminuit>=2
- ipython
- ipywidgets
- joblib
Expand Down
2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@
"bokeh~=1.0",
"eventio>=1.5.0,<2.0.0a0",
"h5py",
"iminuit~=1.3",
"iminuit>=2",
"joblib",
"matplotlib~=3.0",
"numba>=0.43",
Expand Down

0 comments on commit 524e621

Please sign in to comment.