Plot gamma spec from packets

artis-mcrt · Dec 16, 2024 · eb0776e · eb0776e
1 parent 3d9bc4b
commit eb0776e
Show file tree

Hide file tree

Showing 2 changed files with 6 additions and 3 deletions.
diff --git a/artistools/spectra/plotspectra.py b/artistools/spectra/plotspectra.py
@@ -452,6 +452,7 @@ def plot_artis_spectrum(
                 c = 2.99792458e18  # speed of light [angstroms/s]
                 dfspectrum = dfspectrum.with_columns(x=h * c / pl.col("lambda_angstroms") / 1000.0, y=pl.col("f_nu"))
             dfspectrum = dfspectrum.sort("x", maintain_order=True)
+
             axis.plot(
                 dfspectrum["x"], dfspectrum["y"], label=linelabel_withdirbin if axindex == 0 else None, **plotkwargs
             )

diff --git a/artistools/spectra/spectra.py b/artistools/spectra/spectra.py
@@ -59,10 +59,11 @@ def get_exspec_bins(
     mnubins: int | None = None,
     nu_min_r: float | None = None,
     nu_max_r: float | None = None,
+    gamma: bool = False,
 ) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
     """Get the wavelength bins for the emergent spectrum."""
     if modelpath is not None:
-        dfspec = read_spec(modelpath)
+        dfspec = read_spec(modelpath, gamma=gamma)
         if mnubins is None:
             mnubins = dfspec.height
 
@@ -176,7 +177,7 @@ def get_from_packets(
         lambda_bin_edges = np.arange(lambda_min, lambda_max + delta_lambda, delta_lambda)
         lambda_bin_centres = 0.5 * (lambda_bin_edges[:-1] + lambda_bin_edges[1:])  # bin centres
     else:
-        lambda_bin_edges, lambda_bin_centres, delta_lambda = get_exspec_bins(modelpath=modelpath)
+        lambda_bin_edges, lambda_bin_centres, delta_lambda = get_exspec_bins(modelpath=modelpath, gamma=gamma)
         lambda_min = lambda_bin_centres[0]
         lambda_max = lambda_bin_centres[-1]
 
@@ -340,7 +341,8 @@ def get_from_packets(
             "lambda_angstroms",
             pl.col(f"f_lambda_dirbin{dirbin}").alias("f_lambda"),
             pl.col(f"count_dirbin{dirbin}").alias("packetcount"),
-        ])
+            (299792458.0 / (pl.col("lambda_angstroms") * 1e-10)).alias("nu"),
+        ]).with_columns(f_nu=(pl.col("f_lambda") * pl.col("lambda_angstroms") / pl.col("nu")))
         if nprocs_read_dfpackets is None:
             npkts_selected = dfdict[dirbin].get_column("packetcount").sum()
             print(f"    dirbin {dirbin:2d} plots {npkts_selected:.2e} packets")