Commit f843a0c0 authored by Andrea Giannetti's avatar Andrea Giannetti
Browse files

Fixed typo in analytical expressions argument

parent cd03e039

etl/prs/prs_analytical_representations.py

+7 −7
Original line number Diff line number Diff line
@@ -12,7 +12,7 @@ from typing import Union, List 




def approx(x, x_shift, b, exponent, scale):

    arg = 2 * (x - x_shift) / scale

    arg = (x - x_shift) / scale

    clean_arg = np.where(np.abs(arg) > 0.99999, np.sign(arg) * 0.99999, arg)

    return (np.arctanh(clean_arg) + b) ** exponent
@@ -28,11 +28,11 @@ def main(ratios_to_fit: Union[List[str], None] = None): 
    logger.info(f'Data contains {data.shape[0]} rows')



    x0 = {

        '87-86': ((0.55, 6.75, 0.9, 1),),

        '88-87': ((0.59, 6.9, 1.05, 0.96),),

        '88-86': ((0.55, 6.3, 0.83, 1.08),),

        '257-256': ((5, 4.45, 0.75, 8), (4.5, 8.5, 0.79, -7)),

        '381-380': ((2.04, 5.1, 0.6, 2.1), (2.18, 5.7, 0.79, -2.4)),

        '87-86': ((0.55, 6.75, 0.9, 0.5),),

        '88-87': ((0.59, 6.9, 1.05, 0.48),),

        '88-86': ((0.55, 6.3, 0.83, 0.54),),

        '257-256': ((5, 4.45, 0.75, 4), (4.5, 8.5, 0.79, -3.5)),

        '381-380': ((2.04, 5.1, 0.6, 1.05), (2.18, 5.7, 0.79, -1.2)),

    }

    min_ratio = {

        '87-86': 0.08,
@@ -78,7 +78,7 @@ def main(ratios_to_fit: Union[List[str], None] = None): 
                                          kde_max[1][clean_mask],

                                          np.log10(kde_max[0][clean_mask]),

                                          p0=guesses,

                                          bounds=((0.1, 2, 0.3, -9), (7, 10, 1.5, 9)))

                                          bounds=((0.1, 2, 0.3, -4.5), (7, 10, 1.5, 4.5)))

            approx_density = 10 ** approx(x_reg, *best_fit)

            print(best_fit)

            best_fit_params[ratio_string].append({param: float(element) for element, param in zip(best_fit, ['a', 'b', 'nu', 's'])})

etl/prs/prs_make_profiles_figure.py

+4 −3
Original line number Diff line number Diff line
@@ -9,11 +9,11 @@ q = -0.5 
n_0 = 1e5

T_0 = 25



radius = np.linspace(0, 0.9, 100)

radius = np.linspace(0, 0.9, 200)

density = n_0 * (radius / 0.5) ** p

density[0] = 1e8

density = np.where(density > 1e8, 1e8, density)

temperature = T_0 * (radius / 0.5) ** q

temperature[0] = 2000

temperature[0] = temperature[1]



plt.clf()

fig, ax = plt.subplots(layout='constrained')
@@ -26,6 +26,7 @@ ln2 = ax2.plot(radius, temperature, color='red', label='T') 
lns = ln1 + ln2

labs = [l.get_label() for l in lns]

ax.legend(lns, labs, loc=0)

ax2.set_ylim(10, 500)

ax2.semilogy()

ax2.set_ylabel('Temperature [K]')

plt.savefig(os.path.join('prs', 'output', 'density_temperature_profiles.png'))