62 lines
1.8 KiB
Python
62 lines
1.8 KiB
Python
import numpy as np
|
|
|
|
from .thresh import calc_apply_threshold
|
|
from .utils import npmemo
|
|
|
|
|
|
def calc_radial_integration(results, data, pixel_mask):
|
|
do_radial_integration = results.get("do_radial_integration", False)
|
|
if not do_radial_integration:
|
|
return
|
|
|
|
center_x = results["beam_center_x"]
|
|
center_y = results["beam_center_y"]
|
|
|
|
rad, norm = prepare_radial_profile(data.shape, center_x, center_y, pixel_mask)
|
|
|
|
r_min = min(rad)
|
|
r_max = max(rad) + 1
|
|
|
|
data = calc_apply_threshold(results, data, value=np.nan, copy=True)
|
|
|
|
rp = radial_profile(data, rad, norm, pixel_mask)
|
|
|
|
silent_min = results.get("radial_integration_silent_min", None)
|
|
silent_max = results.get("radial_integration_silent_max", None)
|
|
|
|
if silent_min is not None and silent_max is not None:
|
|
# if start > stop, numpy returns an empty array -- better to ensure start < stop by switching them if needed
|
|
silent_min, silent_max = sorted((silent_min, silent_max))
|
|
if silent_min > r_min and silent_max < r_max:
|
|
silent_region = rp[silent_min:silent_max]
|
|
integral_silent_region = np.sum(silent_region)
|
|
rp = rp / integral_silent_region
|
|
results["radint_normalised"] = [silent_min, silent_max]
|
|
|
|
results["radint_I"] = rp[r_min:] #TODO: why not stop at r_max?
|
|
results["radint_q"] = [r_min, r_max]
|
|
|
|
|
|
|
|
@npmemo
|
|
def prepare_radial_profile(shape, x0, y0, keep_pixels):
|
|
y, x = np.indices(shape)
|
|
rad = np.sqrt((x - x0)**2 + (y - y0)**2)
|
|
if keep_pixels is not None:
|
|
rad = rad[keep_pixels]
|
|
rad = rad.astype(int).ravel()
|
|
norm = np.bincount(rad)
|
|
return rad, norm
|
|
|
|
|
|
def radial_profile(data, rad, norm, keep_pixels):
|
|
if keep_pixels is not None:
|
|
data = data[keep_pixels]
|
|
data = data.ravel()
|
|
tbin = np.bincount(rad, data)
|
|
rp = tbin / norm
|
|
return rp
|
|
|
|
|
|
|