How Can I Use Nuker2D? #208
Replies: 2 comments
-
Hi @PatrickHorlaville, thank you for reaching out! There is a Nuker profile available in PetroFit but recent updates in astropy prevents us from using it to fit in its current form. As you know, the Nuker profile explodes when gamma is set to a positive number because you get a 1/0 division at r=0. To deal with this we masked any pixel at r=0 by assigning a If this is critical for your work, you can define a custom model using the astropy custom models. For example you can truncate the image at a radius close to r=0. That is to say for any points sampling at import numpy as np
from astropy.modeling import custom_model
@custom_model
def TruncatedNuker2D(x, y, amplitude=1, r_break=1, x_0=0, y_0=0,
alpha=2, betta=4, gamma=0, epsilon=1):
"""Two dimensional Nuker2D model with smooth transition"""
x_maj = (x - x_0) + (y - y_0)
x_min = -(x - x_0) + (y - y_0)
r = np.sqrt((x_maj) ** 2 + (x_min) ** 2)
# Apply a smooth transition at small radii
r = np.where(r < epsilon, epsilon, r)
return amplitude * (2 ** ((betta - gamma) / alpha)) * (r / r_break) ** (- gamma) * (1 + (r / r_break) ** alpha) ** ((gamma - betta) / alpha) After doing this you can define and fit your model as follows: import petrofit as pf
model_init = TruncatedNuker2D(
amplitude=10,
r_break=10,
x_0=50,
y_0=50,
alpha=2,
betta=4,
gamma=0,
epsilon=1, # 1 pixel truncation
fixed={'epsilon': True}
)
fitted_model, fit_info = pf.fit_model(
image, model_init,
calc_uncertainties=True,
maxiter=10000,
epsilon=1.4901161193847656e-08,
acc=1e-09,
)
pf.plot_fit(fitted_model, image)
plt.show() Notice that |
Beta Was this translation helpful? Give feedback.
-
Here is a version with elliptical isophots: @custom_model
def TruncatedNuker2D(x, y, amplitude=1, r_break=1, x_0=0, y_0=0,
alpha=2, betta=4, gamma=0, theta=0, ellip=0, epsilon=1):
"""Two dimensional Nuker2D model with smooth transition"""
cos_theta = np.cos(theta)
sin_theta = np.sin(theta)
x_maj = np.abs((x - x_0) * cos_theta + (y - y_0) * sin_theta)
x_min = np.abs(-(x - x_0) * sin_theta + (y - y_0) * cos_theta)
r = np.sqrt((x_maj) ** 2 + (x_min/(1 - ellip)) ** 2)
# Apply a smooth transition at small radii
r = np.where(r < epsilon, epsilon, r)
return amplitude * (2 ** ((betta - gamma) / alpha)) * (r / r_break) ** (- gamma) * (1 + (r / r_break) ** alpha) ** ((gamma - betta) / alpha)
You can also restrict the bounds of the model by including a center_slack = 3 # pixels
model_init = TruncatedNuker2D(
amplitude=10,
r_break=10,
x_0=50,
y_0=50,
alpha=2,
betta=4,
gamma=0,
epsilon=1, # 1 pixel truncation
fixed={'epsilon': True}
bounds = {
'x_0': (x_0 - center_slack/2, x_0 + center_slack/2),
'y_0': (y_0 - center_slack/2, y_0 + center_slack/2),
'gamma': (0, None ),
'alpha': (1e-4, None ),
'betta': (1e-4, None ),
'ellip': (0, 0.99),
'theta': (-6.283185307179586, 6.283185307179586),
'r_break': (epsilon, None)
}
) |
Beta Was this translation helpful? Give feedback.
-
Hi all,
I am trying to use PetroFit to model the light profile of a galaxy using the
Nuker2D
function. I have followed the tutorial on the website (the one under 'Image Fitting') and have successfully fitted aSersic2D
profile to my galaxy, but how can I now fit theNuker2D
profile instead?I noticed that the function is present in the models.py file, however if I try to import it with
models.Nuker2D
, an error pops up indicating that the module doesn't have aNuker2D
attribute. I noticed that the Sérsic function was particular in that it had a class for itself,GenSersic2D
, butNuker2D
didn't. Do I need to write a separate class for the Nuker2D profile? If so, I am not sure of how to do so. If anyone can provide any help, that would be greatly appreciated.Thank you!
Beta Was this translation helpful? Give feedback.
All reactions