Hi,

I'm stumping with this not-so-descriptive error when trying to minimize with minuit on jetset 1.1.2

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
<ipython-input-322-f38f62ed22e3> in <module>
      8                                      1E29,
      9                                      fitname='SSC-best-fit-minuit',
---> 10                                      repeat=3)

~/Software/miniconda/envs/astro/lib/python3.7/site-packages/jetset-1.1.2-py3.7-macosx-10.7-x86_64.egg/jetset/minimizer.py in fit(self, fit_Model, sed_data, nu_fit_start, nu_fit_stop, fitname, fit_workplace, loglog, silent, get_conf_int, max_ev, use_facke_err, use_UL, skip_minimizer, repeat)
    393                     if silent is False:
    394                         print('fit run:',i)
--> 395                 self.minimizer.fit(self,max_ev=max_ev,silent=silent)
    396 
    397                 self.pout = self.minimizer.pout

~/Software/miniconda/envs/astro/lib/python3.7/site-packages/jetset-1.1.2-py3.7-macosx-10.7-x86_64.egg/jetset/minimizer.py in fit(self, model, max_ev, use_UL, silent)
    485         self.molde=model
    486         self.silent=silent
--> 487         self._fit(max_ev)
    488         self._fit_stats()
    489         self._set_fit_errors()

~/Software/miniconda/envs/astro/lib/python3.7/site-packages/jetset-1.1.2-py3.7-macosx-10.7-x86_64.egg/jetset/minimizer.py in _fit(self, max_ev)
    686     def _fit(self,max_ev=None):
    687         bounds = [(par.fit_range_min, par.fit_range_max) for par in self.model.fit_par_free]
--> 688         self._set_minuit_func(self.model.pinit, bounds)
    689         #print('=>,')
    690         if max_ev is None or max_ev==0:

~/Software/miniconda/envs/astro/lib/python3.7/site-packages/jetset-1.1.2-py3.7-macosx-10.7-x86_64.egg/jetset/minimizer.py in _set_minuit_func(self, p_init, bounds, p_error)
    732             pedantic=False,
    733             errordef=1,
--> 734             **kwdarg)
    735 
    736     def chisq_func(self, *p):

~/Software/miniconda/envs/astro/lib/python3.7/site-packages/iminuit/minuit.py in __init__(self, fcn, grad, name, *args, **kwds)
    583         )
    584 
--> 585         self._init_state = _make_init_state(self._pos2var, args, kwds)
    586         self._values = mutil.ValueView(self)
    587         self._errors = mutil.ErrorView(self)

~/Software/miniconda/envs/astro/lib/python3.7/site-packages/iminuit/minuit.py in _make_init_state(pos2var, args, kwds)
   1537             if kw not in pos2var:
   1538                 raise RuntimeError(
-> 1539                     f"{kw} is not one of the parameters [{' '.join(pos2var)}]"
   1540                 )
   1541         nargs = len(kwds)

RuntimeError: forced_parameters is not one of the parameters []

Reading here sncosmo/sncosmo#291 it suggest a change in iminuit api. I'm using iminuit 2.4 .

Hi @andreatramacere,

'Hadronic pp jez model' page on Jetset-v1.2.2 doc page mentions that IC cooling was switch off to speed up the process. And, testing galactic_expanding_shell branch, it seems that IC cooling is still off, is it correct?

I was going through LHASSO J2108+5157 paper on MWL study, where jestet-1.2.2 is used and mentions about taking into account IC. So, just wondering if its based on some other development or I am doing some goof-up!

Thank a lot :-)

I have done a couple of fits:

I want to make the plot a bit nicer, so I explore each model component and plot it with matplotlib. I get:

I really don't understand this. I have all these Jet instances stored in a dict and I'm doing exactly the same eval() thing for all of them, and using the same plotting code. It is not a bug in the code because the weird output happens only sometimes, other times it is plot correctly. I have been dealing with this issue for months and I have found no logical pattern, maybe you can bring some light to it ??

PS. I also have other issues with jetset, like how the fit is done (e.g. why if I set a maximum of calls to 300 it does whatever it wants, being it 20 calls or 500?, or why if I set the starting parameters to a reasonable solution the minimizer goes crazy exploring a random part of the space parameter that is super far away from the initial parameter values) or why if I store the Jet models then recovering them brings me wrong parameters?. But for now let's concentrate on the plotting thing I presented before :P

Hello Andrea,

I am learning to use jetset by working with the user guide. Most of the commands work well without any problem including MCMC fitting. However, I got an error when trying to generate the external photon field energy density along the jet. The input was:
fig = plt.figure(figsize=(8,6))
ax=fig.subplots(1)
N=50
G=1
R_range=np.logspace(13,25,N)
y=np.zeros((8,N))
my_jet.set_verbosity(0)
my_jet.set_par('R_BLR_in',1E17)
my_jet.set_par('R_BLR_out',1.1E17)
for ID,R in enumerate(R_range):
my_jet.set_par('R_H',val=R)
my_jet.set_external_fields()
my_jet.energetic_report(verbose=False)

it threw the following error:
An exception has occurred, use %tb to see the full traceback.
SystemExit: JetkerneltException (line 111): the jetkernel failed
exception message: enegetic name thisown not understood

I also have a physics question: when computing EC with the BLR, does it take into account the Klein-Nishina effect? Also, is the absorption of gamma-rays due to pair production with the BLR photons taken into account for the EC calculation? I guessed it because Donnea03 paper deals mainly with the gamma-ray absorption. If so, could you please add a small example about how one can plot the BLR opacity for gamma-gamma absorption as a function of photon energy or distance along the jet.

Thanks.

Hi @andreatramacere ,
I am Ayon Mondal, a final year Master's student in Astrophysics at Presidency University, Kolkata, India.
I am using Jetset for modeling the multi-epoch SED of Blazars as a part of my Msc thesis work. When I am using a very well data the code works fine. But my multi-epoch data has many gaps in different wavebands. For this the pre-fit part works fine, but, after that in the constrains parameters from observable section it shows N_res was negative nan. How can I solve this ?

I am attaching the screenshots herewith.

Your response in this matter would be much appreciated. I am sincerely waiting for your reply.

Hi Adrea,

Currently I start to learn the custom emitters distribution. So I follow the steps in the page:
Custom emitters distribution

But, when I try to start, I've a module error (my python file is called ElectronDistribution as you can see):

  File "ElectronDistribution.py", line 4, in <module>
    from jetset.jet_emitters import EmittersDistribution
ModuleNotFoundError: No module named 'jetset.jet_emitters'

So, I don't know of this is a error due a bad installation (I follow the instruction for conda in a linux machine), so I will apreciatted any kind of help from yours.

Hope to hear about from you soon.

My best,

Antonio.

Thanks, @andreatramacere ,

So this is the new ecsv I have created photandseds_3c279.ecsv.

I have made changes to lines 60, 62, 64 of jetset.data_loader.py as follows
From
self._names = ['x', 'dx', 'y', 'dy', 'T_start', 'T_stop', 'UL', 'data_set']

self._dt = ('f8', 'f8', 'f8', 'f8', 'f8', 'f8', 'bool', 'S16')

self._units = [u.Hz, u.Hz,(u.erg / (u.cm ** 2 * u.s)) ,(u.erg / (u.cm ** 2 * u.s)), cds.MJD, cds.MJD, None, None]

To
self._names = ['No', 'Observed_Passband', 'Photometry_Measurement', 'Uncertainty', 'Units', 'Frequency', 'Flux_Density', 'Upper_limit_of_uncertainty', 'Lower_limit_of_uncertainty', 'Upper_limit_of_Flux_Density', 'Lower_limit_of_Flux_Density', 'NED_Uncertainty', 'Refcode', 'Significance', 'Published_frequency', 'Frequency_Mode', 'Coordinates_Targeted', 'Spatial_Mode', 'Qualifiers', 'Comments']

self._dt = ('f8', 'S16', 'f8', 'f8', 'S16', 'f8', 'f8', 'S16', 'S16', 'S16', 'S16', 'f8', 'S16', 'S16', 'S16', 'S16', 'S16', 'S16', 'S16', 'S16')

self._units = [None, None, u.erg/(u.cm2 * u.s), u.erg/(u.cm2 * u.s), None, u.Hz, u.erg/(u.cm2 * u.s * u.Hz), u.erg/(u.cm2 * u.s * u.Hz), u.erg/(u.cm2 * u.s * u.Hz), u.erg/(u.cm2 * u.s * u.Hz), u.erg/(u.cm2 * u.s * u.Hz), u.erg/(u.cm2 * u.s), None, None, u.Hz, None, None, None, None, None]
So far, no complaints when I ran these (above) lines
There is a complaint on line 113 corresponding to the method from_file().

Further proposed changes are on
lines: 122 - 140 _check_frame_and_scale()
change: all previous labels 'x', 'dx', 'y', 'dy', etc...
to: new labels, e.g. 'No', 'Observed_Passband', 'Photometry_Measurement', etc...

lines: 237 - 244 from_asdc()
change: all previous labels 'x', 'dx', 'y', 'dy', etc...
to: new labels, e.g. 'No', 'Observed_Passband', 'Photometry_Measurement', etc...

lines: 447 - 463 _build_data()
change: instances of 'x', 'y', 'dx', 'dy'
to: 'Frequency', 'Flux_Density', 'Uncertainty'

Or maybe instead of replacing the existing headings, I could simply expand the list by including ones in the photandseds_3c279.ecsv file that are not represented by code?

I hope this make sense?

Hi,
Thanks for this very nice tool.
Just a question about the host galaxy template. If I want to upload a template from a file, where should I give it to the model? For example, when I do

mm,best_fit=my_shape.sync_fit(check_host_gal_template=True,
Ep_start=None,
minimizer='lsb',
silent=True,
fit_range=[10. , 21.])

my_shape.IC_fit(fit_range=[23., 29.],minimizer='minuit',silent=True)
p=my_shape.plot_shape_fit()
p.rescale(y_min=-15)

It checks the host galaxy and returns the parameters. But I want to give the template from a file.
Can you please help with this ?

I cannot load a dataset into the ObsData object, not even the default datasets in jetset.test_data_helper.
Here a snippet

from jetset.data_loader import ObsData
from astropy.table import Table
from jetset.test_data_helper import  test_SEDs
 
mySED = test_SEDs[1]                                                        
sed_data = ObsData(data_file=mySED) 

gives me the error

Traceback (most recent call last):
  File "test_loader.py", line 6, in <module>
    sed_data = ObsData(data_file=mySED)
  File "/home/cosimo/software/miniconda3/lib/python3.7/site-packages/jetset-1.0.2-py3.7-linux-x86_64.egg/jetset/data_loader.py", line 265, in __init__
    raise RuntimeError('table is not an astropy Table')
RuntimeError: table is not an astropy Table

looking into the code
https://github.com/andreatramacere/jetset/blob/master/jetset/data_loader.py#L223
ObsData has a data_table argument in the __init__ than then is checked wheter it is a Table instance
https://github.com/andreatramacere/jetset/blob/master/jetset/data_loader.py#L262

But even generating a dummy astropy table

table_string = """# %ECSV 0.9
# ---
# datatype:
# - {name: x, unit: Hz, datatype: float64}
# - {name: y, unit: erg / (cm2 s), datatype: float64}
# - {name: dy, unit: erg / (cm2 s), datatype: float64}
# schema: astropy-2.0
x y dy
6.204e-05 6.555e-13 1.74e-13
0.000153032 1.4282e-12 4.255e-13
0.00039292 2.052e-12 1.14e-13
0.000355696 2.7004e-12 7.31e-13
0.000947144 4.0304e-12 9.389e-13
0.570768 1.00878e-11 3.7122e-12
0.756888 7.7958e-12 4.392e-12"""

t = Table.read(table_string, format='ascii')
print(t)  

sed_data = ObsData(data_table=t, data_scale="lin-lin", z=1) 

I cannot load it in ObsData

Traceback (most recent call last):
  File "test_loader.py", line 28, in <module>
    sed_data = ObsData(data_table=t, data_scale="lin-lin", z=1) 
  File "/home/cosimo/software/miniconda3/lib/python3.7/site-packages/jetset-1.0.2-py3.7-linux-x86_64.egg/jetset/data_loader.py", line 284, in __init__
    self._build_data(dupl_filter=dupl_filter)
  File "/home/cosimo/software/miniconda3/lib/python3.7/site-packages/jetset-1.0.2-py3.7-linux-x86_64.egg/jetset/data_loader.py", line 447, in _build_data
    self._set_data_frame_and_scale()
  File "/home/cosimo/software/miniconda3/lib/python3.7/site-packages/jetset-1.0.2-py3.7-linux-x86_64.egg/jetset/data_loader.py", line 510, in _set_data_frame_and_scale
    self.data['nu_data_log']+= np.log10(nu_conv_factor)
UnboundLocalError: local variable 'nu_conv_factor' referenced before assignment

Any idea on how to solve it and why not even default dataset work?

Hi,
I am Krishna, a phd scholar from Manipal, Karnataka, India, working in research area of AGN.
I wanted to know that how can I save the fitted SED values (theoretical model lines) into a txt file in JetSet?
so that one can plot else where. Your response in this matter would be much
appreciated. I am sincerely waiting for your reply.

Hello @andreatramacere,

how do I fix the distance of the blob from the target photon field?
For example if I want to consider the external Compton on the dusty torus I get, from the show_pars():

model parameters:
 Name             | Type                 | Units            | value         | phys. boundaries              | log
-------------------------------------------------------------------------------------------------------------------
 B                | magnetic_field       | G                | +1.000000e+00 | [+0.000000e+00,No           ] | False 
 N                | electron_density     | cm^-3            | +1.302546e+03 | [+0.000000e+00,No           ] | False 
 R                | region_size          | cm               | +1.600000e+01 | [+0.000000e+00,+3.000000e+01] | True 
 R_DT             | DT                   | cm               | +5.000000e+18 | [+0.000000e+00,No           ] | False 
 T_DT             | DT                   | K                | +1.000000e+02 | [+0.000000e+00,No           ] | False 
 beam_obj         | beaming              |                  | +1.000000e+01 | [+1.000000e+00,No           ] | False 
 gmax             | high-energy-cut-off  | Lorentz-factor   | +1.000000e+05 | [+1.000000e+00,+1.000000e+15] | False 
 gmin             | low-energy-cut-off   | Lorentz-factor   | +1.000000e+02 | [+1.000000e+00,+1.000000e+05] | False 
 p                | HE_spectral_slope    |                  | +2.800000e+00 | [-1.000000e+01,+1.000000e+01] | False 
 tau_DT           | DT                   |                  | +1.000000e-01 | [+0.000000e+00,+1.000000e+00] | False 
 z_cosm           | redshift             |                  | +6.959965e-02 | [+0.000000e+00,No           ] | False 
-------------------------------------------------------------------------------------------------------------------

How do I evaluate the spectrum from different distances of the blob from the torus?
The Broad Line Region seems to have this R_H

 R_H              | Disk                 | cm               | +1.000000e+17 | [+0.000000e+00,No           ] | False

But I do not find the equivalent for the Torus.

Cheers!

Dear Andrea,
I have installed JetSet from anaconda on my Linux machine (Ubuntu 18 and python 3.8) and tested it successfully following the installation guide. I am now going trough the tutorials and I get stuck in the "QuickStart_1" notebook because of the following error:

ImportError Traceback (most recent call last)
in
----> 1 from jetset.jet_model import Jet
2 my_jet=Jet(name='my_jet',electron_distribution='lppl',)

~/Software/jetset-stable/jetset/jet_model.py in
15 import warnings
16
---> 17 from .jet_spectral_components import JetSpecComponent, SpecCompList
18
19 on_rtd = os.environ.get('READTHEDOCS', None) == 'True'

~/Software/jetset-stable/jetset/jet_spectral_components.py in
24 from .mock import jetkernel as BlazarSED
25 else:
---> 26 from jetkernel import jetkernel as BlazarSED
27
28 from . import spectral_shapes

~/Software/jetset-stable/jetkernel/jetkernel.py in
11 # Import the low-level C/C++ module
12 if package or "." in name:
---> 13 from . import _jetkernel
14 else:
15 import _jetkernel

ImportError: cannot import name '_jetkernel' from 'jetkernel' (/home/marina/Software/jetset-stable/jetkernel/init.py)

Do you maybe have a suggestion?
Thanks a lot in advance!
Marina

Hello @andreatramacere,

I would like to generate an SED with a fixed grid in nu, now if I try to set a value lower than 200 it is just ignored and
an SED of 200 points is returned.
Here a snippet:

from jetset.jet_model import Jet

pwl_jet = Jet(name="", electron_distribution="pl")
pwl_jet.set_nu_grid_size(50)
pwl_jet.show_model()
pwl_jet.eval()
pwl_synch_sed = pwl_jet.get_spectral_component_by_name("Sync").get_SED_points()
# check how many SED points the we get
print(len(pwl_synch_sed[0]))

I get

200

despite in the show_model overview 50 is correctly reported

SED info:
 nu grid size :50
 nu mix (Hz): 1.000000e+06
 nu max (Hz): 1.000000e+30

Is 200 a hardcoded minimum you set?
Thank you!

Hi, I just installed jetset with anaconda without any issues. But when I run the "pytest" command I get the following error:-

pytest --disable-warnings --pyargs jetset.tests.test_users::test_short
============================= test session starts ==============================
platform linux -- Python 3.8.12, pytest-6.2.5, py-1.11.0, pluggy-1.0.0
rootdir: /home/username
collected 0 items

============================ no tests ran in 0.00s =============================
ERROR: module or package not found: jetset.tests.test_users::test_short (missing init.py?)

I did have jetset environment activated before I ran the test.

Secondly, I do not understand the instruction about "run all the examples outside of the installation dir" when installing with anaconda. As per the previous instructions, we created a separate environment called jetset and installed jetset over there, but we did not create any installation directory in the sense we normally would when installing by downloading the tar file->unzipping-> . . .? So, exactly what is the directory we are supposed to exit before running the examples? I am sorry if I am missing something obvious but could you please guide me about this.

P.S.- Best wishes for a merry Christmas and boxing day to this amazing community !

Dear @andreatramacere,

I was making some comparison plots between agnpy, jetset (version 1.1.2) and some reference data I had.
The SSC of jetset always shows a cutoff at high frequencies, compared to agnpy and to the reference SED.

Can you clarify to what it is due and how can I obtain a smoother SSC SED?
Here the snippet of code I use to compute that SED - parameters are set accordingly to the reference SED:

# generate a synchrotron + SSC SED to be confronted with the one produced by
# agnpy and Figure 7.4 of Dermer 2009
from jetset.jet_model import Jet
import numpy as np
import astropy.units as u
import matplotlib.pyplot as plt

# jet with power-law electron distribution
pwl_jet = Jet(
    name="", electron_distribution="pl", electron_distribution_log_values=False
)
# set parameters according to Figure 7.4 of Dermer 2009
pwl_jet.set_par("N", 1298.13238394)
pwl_jet.set_par("p", 2.8)
pwl_jet.set_par("gmin", 1e2)
pwl_jet.set_par("gmax", 1e5)
pwl_jet.set_par("B", 1)
pwl_jet.set_par("R", 1e16)
pwl_jet.set_par("beam_obj", 10)
pwl_jet.set_par("z_cosm", 0.07)
# remove SSA
pwl_jet.spectral_components.Sync.state = "on"
# SSC emission
pwl_jet.set_nu_grid(1e14, 1e26, 100)
pwl_jet.show_model()
pwl_jet.eval()

ssc_nu = pwl_jet.spectral_components.SSC.SED.nu
ssc_sed = pwl_jet.spectral_components.SSC.SED.nuFnu
plt.loglog(ssc_nu, ssc_sed)
plt.ylim([1e-20, 1e-9])
plt.show()

is there some cutoff you are imposing in the numerical calculation?

As always, thanks a lot!

P.S. maybe "cutoff" is not the proper term, I just mean that the SED flux goes abruptly to zero after a given frequency.

Hi @andreatramacere,

Very glad to see the JetSet evolving so fast and now also having preliminary implementation for Galactic objects! 👏🏻

Isn't that the very first line on the documentation page should read as 1) GalacticBeamed to hadle galatcit objects that are not beamed 2) GalacticUnbeamed to hadle galatcit objects that are not beamed ?

Thanks!

Dear Andrea,

I am using jet set for the broadband spectral modelling of a Blazer source with Synchrotron + IC, but I am getting the following error:

*** Log-Polynomial fitting of the IC component ***
---> fit range: [23.0, 29.0]
---> LogCubic fit
---> LogCubic fit failed forced_parameters is not one of the parameters []
---> try LogParabola

RuntimeError Traceback (most recent call last)
~/anaconda3/envs/jetset/lib/python3.7/site-packages/jetset-1.1.2-py3.7-linux-x86_64.egg/jetset/sed_shaper.py in IC_fit(self, fit_range, use_log_par, Ep_start, minimizer, silent)
878 try:
--> 879 mm,best_fit=fit_SED(self.IC_fit_model,self.sed_data,10.**fit_range[0],10.**fit_range[1],loglog=True,silent=True,fitname='IC-shape-fit',minimizer=minimizer)
880
..................
..................

~/anaconda3/envs/jetset/lib/python3.7/site-packages/iminuit/minuit.py in _make_init_state(pos2var, args, kwds)
1966 if kw not in pos2var:
1967 raise RuntimeError(
-> 1968 f"{kw} is not one of the parameters [{' '.join(pos2var)}]"
1969 )
1970 nargs = len(kwds)

RuntimeError: forced_parameters is not one of the parameters []

During the piecewise fitting procedure, the Synchrotron segment is not showing any problem, but in IC fit I am getting this issue. What could be the reason and what could be the resolution? Could you kindly suggest?

After installing jetset in the new conda environment by following the installation instructions, I am facing the following issue while importing the Jet class.

MacOS version: 12.7.1
NumPy version: 1.21.6
Numba version: 0.58.0
Jetset version: 1.2.2

---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
Cell In[2], line 1
----> 1 from jetset.jet_model import Jet
      2 my_jet = Jet(name='test', electron_distribution='lppl',)

File ~/anaconda3/envs/jetset/lib/python3.9/site-packages/jetset/jet_model.py:25
     23 from .jet_paramters import *
     24 from .jet_emitters import *
---> 25 from .jet_emitters_factory import EmittersFactory
     26 from .jet_tools import *
     27 from .mathkernel_helper import bessel_table_file_path

File ~/anaconda3/envs/jetset/lib/python3.9/site-packages/jetset/jet_emitters_factory.py:4
      1 __author__ = "Andrea Tramacere"
      3 import numpy as np
----> 4 import numba as nb
      5 import pprint
      6 from jetset.jet_emitters import EmittersDistribution, InjEmittersDistribution

File ~/anaconda3/envs/jetset/lib/python3.9/site-packages/numba/__init__.py:55
     50             msg = ("Numba requires SciPy version 1.0 or greater. Got SciPy "
     51                    f"{scipy.__version__}.")
     52             raise ImportError(msg)
---> 55 _ensure_critical_deps()
     56 # END DO NOT MOVE
     57 # ---------------------- WARNING WARNING WARNING ----------------------------
     60 from ._version import get_versions

File ~/anaconda3/envs/jetset/lib/python3.9/site-packages/numba/__init__.py:40, in _ensure_critical_deps()
     37 if numpy_version < (1, 22):
     38     msg = (f"Numba needs NumPy 1.22 or greater. Got NumPy "
     39            f"{numpy_version[0]}.{numpy_version[1]}.")
---> 40     raise ImportError(msg)
     41 elif numpy_version > (1, 25):
     42     raise ImportError("Numba needs NumPy 1.25 or less")

ImportError: Numba needs NumPy 1.22 or greater. Got NumPy 1.21.

I tried installing a different version of NumPy. But the pip pointed out the following dependency conflicts:

jetset 1.2.2 requires numpy<1.22,>=1.18, but you have numpy 1.26.2 which is incompatible.
numba 0.58.0 requires numpy<1.26,>=1.21, but you have numpy 1.26.2 which is incompatible.

This suggests numba 0.58.0 would support numpy 1.21.6, so the above ImportError is puzzling. Request any thoughts or suggestions regarding this issue.