Some bug fixed and naming changes

analysis/dip_path_plot.py

@@ -1,7 +1,7 @@
 import numpy as np
 from matplotlib import gridspec
 from charged_shells.rotational_path import PairRotationalPath, PathEnergyPlot
-from charged_shells import expansion, interactions, charge_distributions
+from charged_shells import interactions, charge_distributions
 from charged_shells.parameters import ModelParams
 from config import *
 from plot_settings import *
@@ -134,7 +134,7 @@ def model_comparison(save_as=None, save_data=False):
     # pp_energy = interactions.charged_shell_energy(ex1new, ex2new, params)
     # print(f'PP energy: {pp_energy}')
 
-    # Emanuele data
+    # ICi data
     em_data = np.load(ICI_DATA_PATH.joinpath("FIG_5_JANUS").joinpath("FIG_5_JANUS_A").joinpath("pathway.npz"))['arr_0']
     # em_data = np.load(ICI_DATA_PATH.joinpath("FIG_7").joinpath("pathway.npz"))['arr_0']
     # em_data_path = (ICI_DATA_PATH.joinpath("FIG_8").joinpath("FIXEDCHARGE")
@@ -427,7 +427,7 @@ if __name__ == '__main__':
     #                   )
 
     # model_comparison(
-    #                  # save_as=FIGURES_PATH.joinpath("Emanuele_data").joinpath('IC_CS_janus_path.pdf')
+    #                  # save_as=FIGURES_PATH.joinpath("ICi_data").joinpath('IC_CS_janus_path.pdf')
     #                  )
 
     # combined_kappaR_dependence(kappaR=[1, 3, 10], abar=0.5,

analysis/energy_gap.py

@@ -289,7 +289,7 @@ def main():
     # IC_gap_abar()
 
     IC_gap_charge_at_abar(0.3, which_change='changezc', eta_min=-2, eta_max=2,
-                          save_as=FIGURES_PATH.joinpath('Emanuele_data').joinpath('IC_full_amplitude_charge_abar03.png')
+                          save_as=FIGURES_PATH.joinpath('ICi_data').joinpath('IC_full_amplitude_charge_abar03.png')
                           )
 
 

analysis/higher_multipole_matching.py

@@ -33,7 +33,7 @@ def get_kappaR_higher_multipole_dicts(peak: Peak, params: ModelParams, emanuele_
             idx, = np.nonzero(inverse == i)
             relevant_kR_indices = (emanuele_data[idx, 2] <= max_kappaR) * (min_kappa <= emanuele_data[idx, 2])
             kR = emanuele_data[idx, 2][relevant_kR_indices]
-            en = emanuele_data[idx, peak.emanuele_data_column][relevant_kR_indices]
+            en = emanuele_data[idx, peak.ici_data_column][relevant_kR_indices]
             sort = np.argsort(kR)
             if peak.log_y:
                 data_dict_ic[d] = np.stack((kR[sort], np.abs(en)[sort])).T
@@ -122,7 +122,7 @@ if __name__ == '__main__':
 
     dist = [3.]
     IC_peak_energy_plot(abar=0.5, dist=dist, which='pp', min_kappaR=0.1, max_kappaR=50,
-                        save_as=FIGURES_PATH.joinpath("Emanuele_data").joinpath("peak_pp_kappaR_higher_correction_abar05_rho30.png"),
+                        save_as=FIGURES_PATH.joinpath("ICi_data").joinpath("peak_pp_kappaR_higher_correction_abar05_rho30.png"),
                         save_data=False,
                         quad_correction=False,
                         log_y=False

analysis/patch_shape_comparison.py

@@ -232,7 +232,7 @@ def ic_cs_comparison2(save_data=False):
     # ax.plot(theta, potential_ic.T, label='IC', ls=':', linewidth=2, marker='o', markevery=50, mfc='none')
     ax.tick_params(which='both', direction='in', top=True, right=True, labelsize=11)
     ax.set_xlabel(r'$\theta$', fontsize=11)
-    ax.set_ylabel(r'$\phi [mV]$', fontsize=11)
+    ax.set_ylabel(r'$\psi [mV]$', fontsize=11)
     custom_ticks = [0, np.pi / 4, np.pi / 2, 3 * np.pi / 4, np.pi]
     custom_labels = ['$0$', r'$\pi/4$', r'$\pi/2$', r'$3\pi/4$', r'$\pi$']
     plt.axhline(y=0, color='black', linestyle=':')

analysis/peak_heigth.py

@@ -19,7 +19,7 @@ class Peak:
     y_label: str
     ex1: Expansion
     ex2: Expansion
-    emanuele_data_column: int
+    ici_data_column: int
     log_y: bool
     kappaR_axis_in_expansion: int = None
 
@@ -27,7 +27,7 @@ class Peak:
 class PeakEP(Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 4
+        self.ici_data_column = 4
         self.y_label = r'$|V_\mathrm{EP}|$' if log_y else r'$V_\mathrm{EP}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone()
@@ -39,7 +39,7 @@ class PeakEP(Peak):
 class PeakPP(Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 3
+        self.ici_data_column = 3
         self.y_label = r'$V_\mathrm{PP}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone()
@@ -50,7 +50,7 @@ class PeakPP(Peak):
 class PeakSEP(Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 5
+        self.ici_data_column = 5
         self.y_label = r'$|V_\mathrm{sEP}|$' if log_y else r'$V_\mathrm{sEP}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone()
@@ -75,7 +75,7 @@ def get_charge_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Arra
         if ab in np.around(abar_cs, 5):
             idx, = np.nonzero(inverse == i)
             charge = emanuele_data[idx, 0]
-            en = emanuele_data[idx, peak.emanuele_data_column]
+            en = emanuele_data[idx, peak.ici_data_column]
             sort = np.argsort(charge)
             if peak.log_y:
                 data_dict_ic[ab] = np.stack((charge[sort] / 280 + 2, np.abs(en)[sort])).T
@@ -88,9 +88,9 @@ def get_charge_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Arra
     return data_dict_ic, data_dict_cs
 
 
-def get_kappaR_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Array,
+def get_kappaR_energy_dicts(peak: Peak, params: ModelParams, ici_data: Array,
                             kappaR: Array, abar_cs: list, max_kappaR: float = 50, only_loaded: bool = False):
-    abar_ic, inverse = np.unique(emanuele_data[:, 1], return_inverse=True)
+    abar_ic, inverse = np.unique(ici_data[:, 1], return_inverse=True)
 
     energy_fn = mapping.parameter_map_two_expansions(partial(interactions.charged_shell_energy, dist=2),
                                                      peak.kappaR_axis_in_expansion)
@@ -106,8 +106,8 @@ def get_kappaR_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Arra
         ab = np.around(abar_ic[i], 5)
         if ab in np.around(abar_cs, 5):
             idx, = np.nonzero(inverse == i)
-            kR = emanuele_data[idx, 2][emanuele_data[idx, 2] <= max_kappaR]
-            en = emanuele_data[idx, peak.emanuele_data_column][emanuele_data[idx, 1] <= max_kappaR]
+            kR = ici_data[idx, 2][ici_data[idx, 2] <= max_kappaR]
+            en = ici_data[idx, peak.ici_data_column][ici_data[idx, 1] <= max_kappaR]
             sort = np.argsort(kR)
             if peak.log_y:
                 data_dict_ic[ab] = np.stack((kR[sort], np.abs(en)[sort])).T
@@ -120,9 +120,9 @@ def get_kappaR_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Arra
     return data_dict_ic, data_dict_cs
 
 
-def get_abar_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Array, a_bar: Array, kR_cs: list, max_abar: float = 0.8):
+def get_abar_energy_dicts(peak: Peak, params: ModelParams, ici_data: Array, a_bar: Array, kR_cs: list, max_abar: float = 0.8):
 
-    kR_ic, inverse = np.unique(emanuele_data[:, 2], return_inverse=True)
+    kR_ic, inverse = np.unique(ici_data[:, 2], return_inverse=True)
 
     energy_fn = mapping.parameter_map_two_expansions(partial(interactions.charged_shell_energy, dist=2),
                                                      peak.kappaR_axis_in_expansion)
@@ -135,8 +135,8 @@ def get_abar_energy_dicts(peak: Peak, params: ModelParams, emanuele_data: Array,
         kr = np.around(kR_ic[i], 5)
         if kr in np.around(kR_cs, 5):
             idx, = np.nonzero(inverse == i)
-            abar = emanuele_data[idx, 1][emanuele_data[idx, 1] <= max_abar]
-            en = emanuele_data[idx, peak.emanuele_data_column][emanuele_data[idx, 1] <= max_abar]
+            abar = ici_data[idx, 1][ici_data[idx, 1] <= max_abar]
+            en = ici_data[idx, peak.ici_data_column][ici_data[idx, 1] <= max_abar]
             sort = np.argsort(abar)
             if peak.log_y:
                 data_dict_ic[kr] = np.stack((abar[sort], np.abs(en)[sort])).T
@@ -482,7 +482,7 @@ def main():
 
     a_bar = [0.2, 0.5, 0.8]
     # IC_peak_energy_plot(a_bar=a_bar, which='pp', min_kappaR=0.01, max_kappaR=50,
-    #                     save_as=FIGURES_PATH.joinpath('Emanuele_data').joinpath('peak_pp_kappaR_higher_correction_abar02.png'),
+    #                     save_as=FIGURES_PATH.joinpath('ICi_data').joinpath('peak_pp_kappaR_higher_correction_abar02.png'),
     #                     save_data=False,
     #                     quad_correction=False,
     #                     log_y=False
@@ -493,7 +493,7 @@ def main():
 
     a_bar = [0.1, 0.2, 0.3]
     # IC_peak_energy_charge_plot(a_bar=a_bar, which='sep',
-    #                            # save_as=FIGURES_PATH.joinpath('Emanuele_data').joinpath('peak_sep_charge.png'),
+    #                            # save_as=FIGURES_PATH.joinpath('ICi_data').joinpath('peak_sep_charge.png'),
     #                            )
     IC_peak_energy_charge_combined_plot(a_bar,
                                         save_as=FIGURES_PATH.joinpath('final_figures').joinpath('peak_combined_charge.png')
@@ -501,7 +501,7 @@ def main():
 
     kappaR = [0.01, 3.02407, 30]
     # IC_peak_energy_abar_plot(kappaR=kappaR, which='sep', min_abar=0.2, max_abar=0.8,
-    #                     save_as=FIGURES_PATH.joinpath('Emanuele_data').joinpath('peak_sep_abar.png'),
+    #                     save_as=FIGURES_PATH.joinpath('ICi_data').joinpath('peak_sep_abar.png'),
     #                     save_data=False,
     #                     quad_correction=False,
     #                     log_y=False

analysis/peak_heigth_janus.py

@@ -15,7 +15,7 @@ Expansion = expansion.Expansion
 class JanusPeakPP(peak_heigth.Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 3
+        self.ici_data_column = 3
         self.y_label = r'$V_\mathrm{PP, 1}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone()
@@ -28,7 +28,7 @@ class JanusPeakPP(peak_heigth.Peak):
 class JanusPeakPPinv(peak_heigth.Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 4
+        self.ici_data_column = 4
         self.y_label = r'$V_\mathrm{PP, 2}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone().inverse_sign(exclude_00=True)
@@ -41,7 +41,7 @@ class JanusPeakPPinv(peak_heigth.Peak):
 class JanusPeakEP(peak_heigth.Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 4
+        self.ici_data_column = 4
         self.y_label = r'$V_{EP, 1}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone()
@@ -54,7 +54,7 @@ class JanusPeakEP(peak_heigth.Peak):
 class JanusPeakEPinv(peak_heigth.Peak):
 
     def __init__(self, ex: Expansion, log_y: bool = False, kappaR_axis_in_expansion: int = None):
-        self.emanuele_data_column = 4
+        self.ici_data_column = 4
         self.y_label = r'$V_{EP, 2}$'
         self.ex1 = ex.clone()
         self.ex2 = ex.clone().inverse_sign(exclude_00=True)
@@ -64,8 +64,8 @@ class JanusPeakEPinv(peak_heigth.Peak):
         self.kappaR_axis_in_expansion = kappaR_axis_in_expansion
 
 
-def get_charge_energy_dicts(peak: peak_heigth.Peak, params: ModelParams, emanuele_data: Array, sigma0: Array, abar_cs: list, sigma_tilde=0.001):
-    abar_ic, inverse, counts = np.unique(emanuele_data[:, 1], return_counts=True, return_inverse=True)
+def get_charge_energy_dicts(peak: peak_heigth.Peak, params: ModelParams, ici_data: Array, sigma0: Array, abar_cs: list, sigma_tilde=0.001):
+    abar_ic, inverse, counts = np.unique(ici_data[:, 1], return_counts=True, return_inverse=True)
     # abar_ic = [0.1, 0.2, 0.3]
 
     energy_fn = mapping.parameter_map_two_expansions(partial(interactions.charged_shell_energy, dist=2),
@@ -79,8 +79,8 @@ def get_charge_energy_dicts(peak: peak_heigth.Peak, params: ModelParams, emanuel
         ab = np.around(abar_ic[i], 5)
         if ab in np.around(abar_cs, 5):
             idx, = np.nonzero(inverse == i)
-            charge = emanuele_data[idx, 0]
-            en = emanuele_data[idx, peak.emanuele_data_column]
+            charge = ici_data[idx, 0]
+            en = ici_data[idx, peak.ici_data_column]
             sort = np.argsort(charge)
             if peak.log_y:
                 data_dict_ic[ab] = np.stack((charge[sort], np.abs(en)[sort])).T

analysis/quad_path_plot.py

@@ -506,13 +506,13 @@ def main():
     #                     )
 
     # IC_kappaR_dependence(
-    #                      save_as=FIGURES_PATH.joinpath("Emanuele_data").joinpath('IC_quadrupole_kappaR_dep.png')
+    #                      save_as=FIGURES_PATH.joinpath("ICi_data").joinpath('IC_quadrupole_kappaR_dep.png')
     #                      )
     #
-    # IC_abar_dependence(save_as=FIGURES_PATH.joinpath("Emanuele_data").
+    # IC_abar_dependence(save_as=FIGURES_PATH.joinpath("ICi_data").
     #                    joinpath('IC_quadrupole_abar_dep.png'))
     #
-    # IC_sigma0_dependence(save_as=FIGURES_PATH.joinpath("Emanuele_data").
+    # IC_sigma0_dependence(save_as=FIGURES_PATH.joinpath("ICi_data").
     #                      joinpath('IC_quadrupole_charge_dep_abar05_kappaR3.png'))
 
     # combined_kappaR_dependence(kappaR=[1, 3, 10], abar=0.5,

analysis/sEE_minimum.py

@@ -214,7 +214,7 @@ def main():
     #                   )
 
     IC_kappaR_dependence(which_kappa_lines=[0.01, 3., 5., 10., 50.], file_suffix="_1",
-                         # save_as=FIGURES_PATH.joinpath("Emanuele_data").joinpath('IC_sEE_min_kappaR_abar.png')
+                         # save_as=FIGURES_PATH.joinpath("ICi_data").joinpath('IC_sEE_min_kappaR_abar.png')
                          )
 
 

charged_shells/charge_distributions.py

@@ -31,7 +31,7 @@ def create_mapped_rotsym_expansion(a_bar: Array | float,
     :param a_bar: distance between the center and off center charges
     :param kappaR: screening parameter
     :param sigma_tilde: magnitude of off-center charges / 4pi R^2
-    :param l_max: maximal ell value for the expansion
+    :param l_array:  array of necessary multipole numbers
     :param sigma0: total (mean) charge density
     """
     if not isinstance(sigma0, Array):

charged_shells/expansion.py

@@ -30,6 +30,7 @@ class Expansion:
             raise InvalidExpansion('Array of l values should be unique and sorted.')
         self.coefs = self.coefs.astype(np.complex128)
         self._starting_coefs = np.copy(self.coefs)
+        self._rotations = Quaternion([1., 0., 0., 0.])
 
     def __getitem__(self, item):
         return Expansion(self.l_array, self.coefs[item])