now bayes works

indica/models/background_fit.py

@@ -8,7 +8,7 @@
 from indica.models.plasma import example_run as example_plasma
 from scipy.interpolate import interp1d
 from indica.readers.available_quantities import AVAILABLE_QUANTITIES
-from indica.workflows import run_tomo_1d
+#from indica.workflows import run_tomo_1d
 
 
 def example_run(
@@ -92,4 +92,4 @@ def Bremsstrahlung(
     data["bremsstrahlung"]=(brem)
     return data, brem
 
-run_tomo_1d.pi(10968)
+#run_tomo_1d.pi(10968)

indica/models/diode_filters.py

@@ -33,6 +33,7 @@ def __init__(
         etendue: float = 1.0,
         calibration: float = 2.0e-5,
         instrument_method="get_diode_filters",
+        channel_mask=slice(18,28)
     ):
         """
         Filtered diode diagnostic measuring Bremsstrahlung
@@ -63,6 +64,7 @@ def __init__(
         self.calibration = calibration
         self.instrument_method = instrument_method
         self.quantities = AVAILABLE_QUANTITIES[self.instrument_method]
+        self.channel_mask = channel_mask
 
     def _build_bckc_dictionary(self):
         self.bckc = {}
@@ -116,9 +118,9 @@ def __call__(
         if self.plasma is not None:
             if t is None:
                 t = self.plasma.time_to_calculate
-            Ne = self.plasma.electron_density.interp(t=t)
-            Te = self.plasma.electron_temperature.interp(t=t)
-            Zeff = self.plasma.zeff.interp(t=t).sum("element")
+            Ne = self.plasma.electron_density.sel(t=t)
+            Te = self.plasma.electron_temperature.sel(t=t)
+            Zeff = self.plasma.zeff.sel(t=t).sum("element")
            # print(self.plasma.zeff)
         else:
             if Ne is None or Te is None or Zeff is None:
@@ -154,11 +156,12 @@ def __call__(
             t=t,
             calc_rho=calc_rho,
         )
+        los_integral=los_integral.where((los_integral.channel>self.channel_mask.start)&(los_integral.channel<self.channel_mask.stop))
 
         self.los_integral = los_integral
 
         self._build_bckc_dictionary()
-
+        #print("test1",self.bckc)
         return self.bckc
 
 
@@ -190,7 +193,8 @@ def example_run(pulse: int = None, plasma=None, plot: bool = False):
     model.set_los_transform(los_transform)
     model.set_plasma(plasma)
     bckc = model()
-    print(bckc)
+
+    print("test2", bckc)
 
     if plot:
         it = int(len(plasma.t) / 2)

indica/models/helike_spectroscopy.py

@@ -427,7 +427,6 @@ def __call__(
             self._moment_analysis()
 
         self._build_bckc_dictionary()
-
         return self.bckc
 
 

indica/models/interferometry.py

@@ -87,7 +87,7 @@ def __call__(
         self.los_integral_ne = los_integral_ne
 
         self._build_bckc_dictionary()
-
+        print("INTERFEROMETRY", self.bckc)
         return self.bckc
 
 

indica/workflows/bayes_workflow.py

@@ -184,18 +184,6 @@ def plot_bayes_result(
             figheader=figheader,
             ylabel="Intensity (W m^-2)",
         )
-
-    key = "pi.spectra"
-    if key in blobs.keys():
-        _plot_0d(
-            blobs,
-            key,
-            diag_data,
-            f"{key.replace('.', '_')}.png",
-            figheader=figheader,
-            ylabel="Intensity (W m^-2)",
-        )
-
     key = "efit.wp"
     if key in blobs.keys():
         _plot_0d(

indica/workflows/example_bayes_opt.py

@@ -35,7 +35,7 @@ def run(
         dt=dt,
         main_ion="h",
         impurities=("c",),  #impurities: tuple = ("c", "ar"), impurity_concentration: tuple = (0.02, 0.001),
-        impurity_concentration=(0.02,),
+        impurity_concentration=(0.2,),
         full_run=False,
         n_rad=10,
     )
@@ -51,47 +51,28 @@ def run(
     }
 
     ST40 = read_st40.ReadST40(pulse) 
-    ST40(["xrcs", "smmh1", "pi"])
-
-    los_transform = ST40.binned_data["smmh1"]["ne"].transform
-    smmh1 = Interferometry(name="smmh1")
-    smmh1.set_los_transform(los_transform)
-    smmh1.plasma = plasma
-    los_transform = ST40.binned_data["xrcs"]["te_kw"].transform
-    xrcs = Helike_spectroscopy(name="xrcs", window_masks=[slice(0.3945, 0.3962)])
-    xrcs.set_los_transform(los_transform)
-    xrcs.plasma = plasma
-
-    data_to_read=ST40.binned_data["pi"]["spectra"].sel(channel=slice(18,28))
-    los_transform = data_to_read.transform
-    data_to_read.transform.set_equilibrium(
-        data_to_read.transform.equilibrium
-    )
-    pi = BremsstrahlungDiode(name="pi")
+    ST40(["pi"])
 
+    data_to_read=ST40.binned_data["pi"]["spectra"]
+    los_transform = data_to_read.transform
+    data_to_read.transform.set_equilibrium(data_to_read.transform.equilibrium)
+    pi = BremsstrahlungDiode(name="pi", channel_mask=slice(18, 28))
     pi.set_los_transform(los_transform)
 
-    from indica.models.plasma import example_run as example_plasma
-    pi.set_plasma(example_plasma(pulse=pulse, impurities=("c",), impurity_concentration=(0.02,)))
-    print(pi())
-
-    pi_data=pi()["brightness"].sel(channel=slice(18,28))
+    pi.plasma = plasma
+    pi_data = pi()["brightness"]
 
+    #from indica.models.background_fit import Bremsstrahlung
+    #pi_data = Bremsstrahlung(pulse)[1]
+    #data_measured = Bremsstrahlung(pulse)[1].sel(channel=channels)
+    #data_modelled = example_run(pulse)[2]["brightness"].sel(channel=channels)
 
 
     flat_data = {}
-    flat_data["smmh1.ne"] = (
-        smmh1().pop("ne").expand_dims(dim={"t": [plasma.time_to_calculate]})
-    )
-
-    flat_data["xrcs.spectra"] = (
-        xrcs().pop("spectra").expand_dims(dim={"t": [plasma.time_to_calculate]})
-    )
     flat_data["pi.brightness"] = (
-        pi_data
+        pi_data.expand_dims(dim={"t": [plasma.time_to_calculate]})
     )
 
-
     priors = {
         "Ne_prof.y0": get_uniform(1e19, 8e19),
         "Ne_prof.y1": get_uniform(1e18, 5e18),
@@ -129,17 +110,15 @@ def run(
     bm = BayesModels(
         plasma=plasma,
         data=flat_data,
-        diagnostic_models=[smmh1, xrcs, pi], 
+        diagnostic_models=[pi],
         quant_to_optimise=[
-            "smmh1.ne",
-            "xrcs.spectra",
             "pi.brightness",
         ],
         priors=priors,
     )
 
     ndim = param_names.__len__()
-    nwalkers = 20
+    nwalkers = 50
     start_points = bm.sample_from_priors(param_names, size=nwalkers)
     move = [(emcee.moves.StretchMove(), 1.0), (emcee.moves.DEMove(), 0.0)]
     sampler = emcee.EnsembleSampler(
@@ -166,9 +145,7 @@ def run(
     }
 
     samples = sampler.get_chain(flat=True)
-
     prior_samples = bm.sample_from_priors(param_names, size=int(1e5))
-
     result = {
             "blobs": blob_dict,
             "diag_data": flat_data,
@@ -203,4 +180,4 @@ def run(
         pathname="./plots/"
     elif platform == "win32":
         pathname="C:\\Users\\Aleksandra.Alieva\\Desktop\\Plots\\New\\"
-    run(10607, params, 10, pathname, burn_in=0)
+    run(10607, params, 200, pathname, burn_in=0)

indica/workflows/run_tomo_1d.py

@@ -19,7 +19,8 @@
 set_plot_rcparams("profiles")
 
 
-def sxrc_xy(
+"""
+    def sxrc_xy(
     pulse: int = 10820,
     tstart: float = 0.02,
     tend: float = 0.11,
@@ -179,6 +180,7 @@ def old_camera(
        # tomo.show_reconstruction()
 
     return input_dict
+    """
 
 def pi(
     pulse,
@@ -250,7 +252,7 @@ def pi(
                     f"{pulse}_{instrument}_{timestr}_{i}",
                     save_fig=save_fig,
                 )
-        """
+
         plt.figure()
         surf = data_measured.T.plot()
         set_axis_sci(plot_object=surf)
@@ -311,17 +313,19 @@ def pi(
             debug=debug,
             has_data=has_data,
         )
-    """
-    #tomo = tomo_1D.SXR_tomography(input_dict, reg_level_guess=reg_level_guess)
-    #tomo()
+
+    tomo = tomo_1D.SXR_tomography(input_dict, reg_level_guess=reg_level_guess)
+    tomo()
 
     if plot:
         plt.ioff()
-        #tomo.show_reconstruction()
+        tomo.show_reconstruction()
         plt.show()
     return input_dict
 
-def fake_data(
+pi(10968)
+
+"""def fake_data(
     pulse: int = 9229,
     plasma=None,
     model=None,
@@ -394,3 +398,4 @@ def fake_data(
         tomo.show_reconstruction()
 
     return plasma
+    """