Skip to content

Commit 49d0e4a

Browse files
committed
Cone-beam: refactored Ellipse class, and reorganized visualizer
Signed-off-by: Nicola VIGANO <nicola.vigano@esrf.fr>
1 parent 40e95b9 commit 49d0e4a

3 files changed

Lines changed: 268 additions & 238 deletions

File tree

corrct/alignment/cone_beam.py

Lines changed: 14 additions & 144 deletions
Original file line numberDiff line numberDiff line change
@@ -283,8 +283,8 @@ def __init__(
283283
self._pre_fit()
284284

285285
def _pre_fit(self, use_least_squares: bool = False) -> None:
286-
ell1_acq = fitting.Ellipse(self.points_ell1)
287-
ell2_acq = fitting.Ellipse(self.points_ell2)
286+
self.ell1_acq = fitting.Ellipse(self.points_ell1, least_squares=use_least_squares)
287+
self.ell2_acq = fitting.Ellipse(self.points_ell2, least_squares=use_least_squares)
288288

289289
if self.points_axis is not None:
290290
# Using measured projected center, whenever available
@@ -293,8 +293,8 @@ def _pre_fit(self, use_least_squares: bool = False) -> None:
293293

294294
self.prj_origin_vu = self.points_axis[:, 1]
295295
else:
296-
self.ell1_prj_center_vu = ell1_acq.fit_prj_center(least_squares=use_least_squares)
297-
self.ell2_prj_center_vu = ell2_acq.fit_prj_center(least_squares=use_least_squares)
296+
self.ell1_prj_center_vu = self.ell1_acq.center_vu
297+
self.ell2_prj_center_vu = self.ell2_acq.center_vu
298298

299299
self.prj_origin_vu = None
300300

@@ -315,30 +315,25 @@ def _pre_fit(self, use_least_squares: bool = False) -> None:
315315
self.points_ell2_rot = self.points_ell2.copy()
316316

317317
# Re-instatiate ellipse class, after rotation
318-
ell1_rot = fitting.Ellipse(self.points_ell1_rot)
319-
ell2_rot = fitting.Ellipse(self.points_ell2_rot)
320-
321-
self.ell1_params = ell1_rot.fit_parameters(least_squares=use_least_squares)
322-
self.ell2_params = ell2_rot.fit_parameters(least_squares=use_least_squares)
318+
self.ell1_rot = fitting.Ellipse(self.points_ell1_rot, least_squares=use_least_squares)
319+
self.ell2_rot = fitting.Ellipse(self.points_ell2_rot, least_squares=use_least_squares)
323320

324321
if self.plot_result:
325322
fig, axs = plt.subplots()
326323
axs.plot(self.points_ell1[1, :], self.points_ell1[0, :], "C0--", label="Ellipse 1 - Acquired")
327324
axs.plot(self.points_ell2[1, :], self.points_ell2[0, :], "C1--", label="Ellipse 2 - Acquired")
328325
axs.plot(self.points_ell1_rot[1, :], self.points_ell1_rot[0, :], "C0", label="Ellipse 1 - Rotated")
329326
axs.plot(self.points_ell2_rot[1, :], self.points_ell2_rot[0, :], "C1", label="Ellipse 2 - Rotated")
330-
ell1_acq_params = ell1_acq.fit_parameters(least_squares=use_least_squares)
331-
ell2_acq_params = ell2_acq.fit_parameters(least_squares=use_least_squares)
332-
axs.plot([ell1_acq_params[-1], ell2_acq_params[-1]], [ell1_acq_params[-2], ell2_acq_params[-2]], "C2--")
333-
axs.plot([self.ell1_params[-1], self.ell2_params[-1]], [self.ell1_params[-2], self.ell2_params[-2]], "C2")
327+
axs.plot([self.ell1_acq.u, self.ell2_acq.u], [self.ell1_acq.v, self.ell2_acq.v], "C2--")
328+
axs.plot([self.ell1_rot.u, self.ell2_rot.u], [self.ell1_rot.v, self.ell2_rot.v], "C2")
334329
if self.points_axis is not None:
335330
axs.scatter(self.points_axis[1], self.points_axis[0], c="C2", marker="*", label="Centers - Acquired")
336331
axs.legend()
337332
axs.grid()
338333
fig.tight_layout()
339334
plt.show(block=False)
340335

341-
self.acq_geom.D = self._fit_distance_det2src(self.ell1_params, self.ell2_params)
336+
self.acq_geom.D = self._fit_distance_det2src(self.ell1_rot, self.ell2_rot)
342337

343338
if self.verbose:
344339
print(f"Fitted detector distance from source (pix): {self.acq_geom.D}")
@@ -359,8 +354,8 @@ def fit(self, r: float, e: float = 1) -> ConeBeamGeometry:
359354
ValueError
360355
In case of flipped ellipses.
361356
"""
362-
b1, a1, c1, v1, u1 = self.ell1_params
363-
b2, a2, c2, v2, u2 = self.ell2_params
357+
b1, a1, c1, v1, u1 = self.ell1_rot.parameters
358+
b2, a2, c2, v2, u2 = self.ell2_rot.parameters
364359

365360
sign_z1 = -1
366361
sign_z2 = sign_z1 * -e
@@ -430,11 +425,9 @@ def get_zeta(bk, ak, ck, D, sign_zk) -> float:
430425
return self.acq_geom
431426

432427
@staticmethod
433-
def _fit_distance_det2src(
434-
ellipse_1: Union[ArrayLike, NDArray], ellipse_2: Union[ArrayLike, NDArray], e: float = 1
435-
) -> float:
436-
b1, a1, c1, v1, _ = np.array(ellipse_1)
437-
b2, a2, c2, v2, _ = np.array(ellipse_2)
428+
def _fit_distance_det2src(ellipse_1: fitting.Ellipse, ellipse_2: fitting.Ellipse, e: float = 1) -> float:
429+
b1, a1, c1, v1, _ = ellipse_1.parameters
430+
b2, a2, c2, v2, _ = ellipse_2.parameters
438431

439432
ecc2 = np.sqrt(b2 - c2**2 / a2)
440433
ecc1 = np.sqrt(b1 - c1**2 / a1)
@@ -522,126 +515,3 @@ def tune_acquisition_geometry(
522515
acq_geom_tuned = acq_geom_tuned.update(par_name, min_par)
523516

524517
return acq_geom_tuned
525-
526-
527-
def cm2inch(dims: Union[ArrayLike, NDArray]) -> tuple[float]:
528-
"""Convert cm into inch.
529-
530-
Parameters
531-
----------
532-
dims : Union[ArrayLike, NDArray]
533-
The dimentions of the object in cm
534-
535-
Returns
536-
-------
537-
tuple[float]
538-
The output dimensions in inch
539-
"""
540-
return tuple(np.array(dims) / 2.54)
541-
542-
543-
class MarkerVisualizer:
544-
"""Plotting class to assess the calibration quality."""
545-
546-
def __init__(
547-
self,
548-
fitted_positions_vu: Union[ArrayLike, NDArray],
549-
imgs: NDArray,
550-
disk: NDArray,
551-
ell_params: Union[ArrayLike, NDArray, None] = None,
552-
) -> None:
553-
self.positions_vu = np.array(fitted_positions_vu)
554-
self.imgs = imgs
555-
self.disk = disk
556-
self.global_lims = False
557-
558-
if ell_params is not None:
559-
ell_params = np.array(ell_params)
560-
self.ell_params = ell_params
561-
562-
self.curr_pos = 0
563-
564-
if self.ell_params is not None:
565-
us = np.sort(self.positions_vu[1, :])
566-
self.v_1, self.v_2 = fitting.Ellipse.predict_v(self.ell_params, us)
567-
568-
self.fig, self.axs = plt.subplots(1, 3, figsize=cm2inch([36, 12])) # , sharex=True, sharey=True
569-
self.axs[2].imshow(self.disk)
570-
self.axs[0].set_xlim(0, self.imgs.shape[-1])
571-
self.axs[0].set_ylim(self.imgs.shape[-3], 0)
572-
self.fig.tight_layout()
573-
self.update()
574-
575-
self.fig.canvas.mpl_connect("key_press_event", self._key_event)
576-
self.fig.canvas.mpl_connect("scroll_event", self._scroll_event)
577-
578-
def update(self) -> None:
579-
self.curr_pos = self.curr_pos % self.imgs.shape[-2]
580-
581-
for img in self.axs[0].get_images():
582-
img.remove()
583-
x_lims = self.axs[0].get_xlim()
584-
y_lims = self.axs[0].get_ylim()
585-
self.axs[0].cla()
586-
self.axs[0].set_xlim(x_lims[0], x_lims[1])
587-
self.axs[0].set_ylim(y_lims[0], y_lims[1])
588-
589-
for img in self.axs[1].get_images():
590-
img.remove()
591-
self.axs[1].cla()
592-
593-
self.axs[0].plot(self.positions_vu[1, :], self.positions_vu[0, :], "bo-", markersize=4)
594-
self.axs[0].scatter(self.positions_vu[1, self.curr_pos], self.positions_vu[0, self.curr_pos], c="r")
595-
596-
if self.ell_params is not None:
597-
us = np.sort(self.positions_vu[1, :])
598-
self.axs[0].plot(us, self.v_1, "g")
599-
self.axs[0].plot(us, self.v_2, "g")
600-
self.axs[0].grid()
601-
602-
if self.global_lims:
603-
vmin = self.imgs.min()
604-
vmax = self.imgs.max()
605-
else:
606-
vmin = self.imgs[:, self.curr_pos, :].min()
607-
vmax = self.imgs[:, self.curr_pos, :].max()
608-
609-
img = self.axs[1].imshow(self.imgs[:, self.curr_pos, :], vmin=vmin, vmax=vmax)
610-
self.axs[1].scatter(self.positions_vu[1, self.curr_pos], self.positions_vu[0, self.curr_pos], c="r")
611-
self.axs[1].set_title(f"Range: [{vmin}, {vmax}]")
612-
# plt.colorbar(im, ax=self.axs[1])
613-
self.fig.canvas.draw()
614-
615-
def _key_event(self, evnt) -> None:
616-
if evnt.key == "right":
617-
self.curr_pos += 1
618-
elif evnt.key == "left":
619-
self.curr_pos -= 1
620-
elif evnt.key == "up":
621-
self.curr_pos += 1
622-
elif evnt.key == "down":
623-
self.curr_pos -= 1
624-
elif evnt.key == "pageup":
625-
self.curr_pos += 10
626-
elif evnt.key == "pagedown":
627-
self.curr_pos -= 10
628-
elif evnt.key == "escape":
629-
plt.close(self.fig)
630-
elif evnt.key == "ctrl+l":
631-
self.global_lims = not self.global_lims
632-
else:
633-
print(evnt.key)
634-
return
635-
636-
self.update()
637-
638-
def _scroll_event(self, evnt) -> None:
639-
if evnt.button == "up":
640-
self.curr_pos += 1
641-
elif evnt.button == "down":
642-
self.curr_pos -= 1
643-
else:
644-
print(evnt.key)
645-
return
646-
647-
self.update()

0 commit comments

Comments
 (0)