This is a very simple and straightforward way of demonstrating that higher energy neutrons experience stronger anisotropy.
We will calculate the magnitude of anisotropy as a one-dimensional parameter as:
$$A_{sym} = \frac{W(180^\circ)}{W(90^\circ)}$$This is the same method as Tony Shin is using in his multiplicity counter work, and the same method that several previous papers have used.
In [1]:
import os
import sys
import matplotlib.pyplot as plt
import numpy as np
import imageio
In [2]:
import pandas as pd
In [3]:
import seaborn as sns
sns.set(style='ticks')
In [4]:
sys.path.append('../scripts/')
In [5]:
import bicorr as bicorr
import bicorr_plot as bicorr_plot
import bicorr_sums as bicorr_sums
import bicorr_math as bicorr_math
In [6]:
%load_ext autoreload
%autoreload 2
I will work with the data from the combined data sets in analysis\Cf072115_to_Cf072215b.
In [7]:
det_df = bicorr.load_det_df('../meas_info/det_df_pairs_angles.csv', remove_fc_neighbors = True)
In [8]:
singles_hist, dt_bin_edges_sh, dict_det_to_index, dict_index_to_det = bicorr.load_singles_hist(filepath='../analysis/Cf072115_to_Cf072215b/datap')
In [11]:
npzfile = np.load('../analysis_no_Ethresh_with_CT\Cf072115_to_Cf072215b\datap/bhp_nn_by_pair_1ns.npz')
bhp_nn_pos = npzfile['bhp_nn_pos']
bhp_nn_neg = npzfile['bhp_nn_neg']
dt_bin_edges = npzfile['dt_bin_edges']
emin, emax will be input parameters. I will vary emin.
In [12]:
num_fissions = 2194651200.00
angle_bin_edges = np.arange(8,190,10)
In [13]:
emin = .62
emax = 12
The code will produce the sum across energy ranges that are rounded down from these values to find the bin edges. Where should I get them back from?
Add optional flag to return energies_real. This way it won't break all my old notebooks. Note that energies_real is [emax_real,emin_real].
In [14]:
singles_df, det_df, by_angle_df, energies_real = bicorr_sums.perform_W_calcs(det_df,
dict_index_to_det, singles_hist, dt_bin_edges_sh,
bhp_nn_pos, bhp_nn_neg, dt_bin_edges,
num_fissions, emin, emax, angle_bin_edges, return_real_energies_flag = True)
print(energies_real)
In [15]:
angle_bin_edges
Out[15]:
In [16]:
series_180 = by_angle_df.loc[np.int(np.digitize(180,angle_bin_edges))-1]
In [17]:
series_90 = by_angle_df.loc[np.int(np.digitize(90,angle_bin_edges))-1]
In [18]:
series_180['W']/series_90['W']
Out[18]:
In [27]:
series_min = by_angle_df.loc[by_angle_df['W'].idxmin()]
In [28]:
num = series_180['W']
num_err = series_180['W_err']
denom = series_90['W']
denom_err = series_90['W_err']
minval = series_min['W']
minval_err = series_min['W_err']
Asym, Asym_err = bicorr_math.prop_err_division(num,num_err,denom,denom_err)
Asym_min, Asym_min_err = bicorr_math.prop_err_division(num,num_err,minval,minval_err)
In [29]:
print(Asym, Asym_err)
print(Asym_min, Asym_min_err)
In [30]:
bicorr_sums.calc_Asym(by_angle_df)
Out[30]:
In [32]:
bicorr_sums.calc_Asym(by_angle_df, min_flag=True)
Out[32]:
Create a dataframe for storing these values, because I like dataframes now.
In [33]:
emins = np.arange(0.5,5,.2)
print(emins)
emax = 12
In [34]:
Asym_df = pd.DataFrame(data = emins, columns = {'emin'})
Asym_df['emax'] = emax
Asym_df['Asym'] = np.nan
Asym_df['Asym_err'] = np.nan
Asym_df['Asym_min'] = np.nan
Asym_df['Asym_min_err'] = np.nan
In [35]:
Asym_df.head()
Out[35]:
In [36]:
for index, row in Asym_df.iterrows():
singles_df, det_df, by_angle_df = bicorr_sums.perform_W_calcs(det_df,
dict_index_to_det, singles_hist, dt_bin_edges_sh,
bhp_nn_pos, bhp_nn_neg, dt_bin_edges,
num_fissions, row['emin'], row['emax'], angle_bin_edges)
Asym, Asym_err, Asym_min, Asym_min_err = bicorr_sums.calc_Asym(by_angle_df,min_flag=True)
Asym_df.loc[index,'Asym'] = Asym
Asym_df.loc[index,'Asym_err'] = Asym_err
Asym_df.loc[index,'Asym_min'] = Asym_min
Asym_df.loc[index,'Asym_min_err'] = Asym_min_err
In [37]:
plt.figure(figsize=(4,3))
plt.errorbar(Asym_df['emin'],Asym_df['Asym'],yerr=Asym_df['Asym_err'],fmt='.',color='k')
plt.xlabel('$E_{min}$ (MeV)')
plt.ylabel('$A_{sym}$')
plt.title('Errors from std(W)')
sns.despine(right=True)
bicorr_plot.save_fig_to_folder('Asym_vs_angle_std_W')
plt.show()
If instead I want to propagate errors from W_err (statistical error):
In [39]:
for index, row in Asym_df.iterrows():
singles_df, det_df, by_angle_df = bicorr_sums.perform_W_calcs(det_df,
dict_index_to_det, singles_hist, dt_bin_edges_sh,
bhp_nn_pos, bhp_nn_neg, dt_bin_edges,
num_fissions, row['emin'], row['emax'], angle_bin_edges)
Asym, Asym_err, Asym_min, Asym_min_err = bicorr_sums.calc_Asym(by_angle_df,std_flag=True,min_flag=True)
Asym_df.loc[index,'Asym'] = Asym
Asym_df.loc[index,'Asym_err'] = Asym_err
Asym_df.loc[index,'Asym_min'] = Asym_min
Asym_df.loc[index,'Asym_min_err'] = Asym_min_err
In [40]:
plt.figure(figsize=(4,3))
plt.errorbar(Asym_df['emin'],Asym_df['Asym'],yerr=Asym_df['Asym_err'],fmt='.',color='k')
plt.xlabel('$E_{min}$ (MeV)')
plt.ylabel('$A_{sym}$')
plt.title('Errors from W_err')
sns.despine(right=True)
bicorr_plot.save_fig_to_folder('Asym_vs_angle_W_err')
plt.show()
In [41]:
plt.figure(figsize=(4,3))
plt.errorbar(Asym_df['emin'],Asym_df['Asym_min'],yerr=Asym_df['Asym_min_err'],fmt='.',color='k')
plt.xlabel('$E_{min}$ (MeV)')
plt.ylabel('$A_{sym}$')
plt.title('Errors from W_err')
sns.despine(right=True)
bicorr_plot.save_fig_to_folder('Asym_min_vs_angle_W_err')
plt.show()
In [42]:
Asym_df = bicorr_sums.calc_Asym_vs_emin(det_df,
dict_index_to_det, singles_hist, dt_bin_edges_sh,
bhp_nn_pos, bhp_nn_neg, dt_bin_edges,
num_fissions, emins, emax, angle_bin_edges,
plot_flag=True, save_flag=False)
Asym_df.head()
Out[42]:
I'm updating this to return emin_real and emax_real, which are slightly different than emin and emax due to discrete binning.
In [32]:
Asym_df = bicorr_sums.calc_Asym_vs_emin(det_df,
dict_index_to_det, singles_hist, dt_bin_edges_sh,
bhp_nn_pos, bhp_nn_neg, dt_bin_edges,
num_fissions, emins, emax, angle_bin_edges,
plot_flag=True, save_flag=False)
Asym_df.head()
Out[32]:
In [13]:
e_bin_edges = np.arange(0.5,5,.2)
In [17]:
Asym_df = pd.DataFrame(data = {'emin':e_bin_edges[:-1],'emax':e_bin_edges[1:]})
Asym_df.head()
Out[17]:
I updated this method in a new function called bicorr_sums.calc_Asym_vs_ebin. The rest of the function is the same as calc_Asym_vs_emin.
In [20]:
Asym_df = bicorr_sums.calc_Asym_vs_ebin(det_df,
dict_index_to_det, singles_hist, dt_bin_edges_sh,
bhp_nn_pos, bhp_nn_neg, dt_bin_edges,
num_fissions, e_bin_edges, angle_bin_edges,
plot_flag=True, save_flag=False)
Asym_df.head()
Out[20]:
In [ ]: